Dynamic mode decomposition of nonequilibrium electron-phonon dynamics: accelerating the first-principles real-time Boltzmann equation

Nonequilibrium dynamics governed by electron-phonon (e-ph) interactions plays a key role in electronic devices and spectroscopies and is central to understanding electronic excitations in materials. The real-time Boltzmann transport equation (rt-BTE) with collision processes computed from first principles can describe the coupled dynamics of electrons and atomic vibrations (phonons). Yet, a bottleneck of these simulations is the calculation of e-ph scattering integrals on dense momentum grids at each time step. Here we show a data-driven approach based on dynamic mode decomposition (DMD) that can accelerate the time propagation of the rt-BTE and identify dominant electronic processes. We apply this approach to two case studies, high-field charge transport and ultrafast excited electron relaxation. In both cases, simulating only a short time window of ~10% of the dynamics suffices to predict the dynamics from initial excitation to steady state using DMD extrapolation. Analysis of the momentum-space modes extracted from DMD sheds light on the microscopic mechanisms governing electron relaxation to steady state or equilibrium. The combination of accuracy and efficiency makes our DMD-based method a valuable tool for investigating ultrafast dynamics in a wide range of materials.

First-principles calculations in the time domain provide a microscopic description of nonequilibrium dynamics in materials.These methods propagate in time quantities characterizing the quantum dynamics, such as the time-dependent electron wave function, density [24], density matrix [25] or Green's function [21], and can also access the time-dependent atomic positions and lattice vibrations [17,18].Different schemes are successful in different regimes.For example, coherent electron dynamics on the attosecond time scale can be modeled effectively using time-dependent DFT [22,[26][27][28][29], but that approach is not ideal for modeling phonon dynamics, which occurs on a picosecond time scale [30].
The real-time Boltzmann transport equation (rt-BTE) has emerged as an effective tool for exploring the coupled electron and phonon dynamics from femtosecond to nanosecond timescales [15,18,31].In the rt-BTE, the time-dependent electron populations are obtained by solving a set of integro-differential equations accounting for e-ph scattering processes on dense momentum grids.Following an initial excitation, the rt-BTE is propagated in time to reach thermal equilibrium or steady state in an external field.This scheme employs a femtosecond time step to capture the e-ph scattering processes.However, evaluating the scattering integral at each time step makes the rt-BTE approach computationally demanding even for materials with a handful of atoms in the unit cell.
Data-driven techniques are increasingly employed in materials modeling, both for accelerating computational workflows and to gain physical insight using learning algorithms [32,33].In particular, dynamic mode decomposition (DMD), which was developed in the last decade to study fluid dynamics, is a valuable tool to linearize dynamical problems and reduce their dimensionality [34,35].In DMD, explicit simulation of a short initial time window allows one to learn the dominant modes governing the dynamics and extrapolate the simulation to future times at low computational cost.Recent work has employed DMD to study electron dynamics described by model Hamiltonians with purely electronic interactions [36,37].Yet, to date DMD has not been applied to more computationally intensive first-principles studies.
In this work, we combine DMD with first-principles calculations of nonequilibrium electron dynamics, using the framework of the rt-BTE in the presence of e-ph collisions and external fields.We show that DMD provides an order-of-magnitude computational speed-up while retaining the full accuracy of the first-principles rt-BTE.In addition, DMD reveals key momentum-space temporal patterns and achieves a significant dimensionality reduction of the nonequilibrium physics.Our results include both high-field transport and transient excited-state dynamics, and are accompanied by a careful characterization of convergence with respect to the size of the sampling window during which DMD learns the dominant modes.Taken together, this work provides the blueprint for combining data-driven methods with first-principles calculations to study nonequilibrium dynamics in real materials.

A. First-principles rt-BTE
We describe the electron distribution using the timedependent populations f nk (t), which quantify the occupation of each electronic state |nk⟩, where k is the electron crystal momentum and n is the band index (from now on we omit the band index to simplify the notation).Starting from an initial distribution at time zero, f k (t = 0), in the rt-BTE the populations evolve according to [38] where I[f k (t)] is the collision integral accounting for eph scattering processes in momentum space [39] and F includes any external fields applied to the system.The rt-BTE simulations use dense momentum grids to accurately describe scattering between electronic states via absorption and emission of phonons.The required grid sizes are typically greater than 100 × 100 × 100 for both electron and phonon momenta.We time-step equation 1 using explicit solvers (Euler or 4th-order Runge-Kutta) or more advanced Strang splitting techniques [31].The collision integral includes a summation over the phonon momentum grid and is evaluated at least once per time step using a parallel algorithm implemented in the Perturbo code [38] (see Methods for details).Although here we focus on the dynamics of electrons interacting with phonons, the rt-BTE formalism has also been extended to study nonequilibrium phonon [17] and exciton dynamics [40].

B. DMD learning and prediction of the dynamics
We employ DMD in combination with rt-BTE simulations.The DMD approach linearizes the dynamics by relating the states of the system at times t and t + ∆t via a time-independent matrix A [41,42].Focusing on the e-ph dynamics, this amounts to advancing the electronic populations at time t using where the populations f k form a vector with size N equal to the number of k-points in the electronic momentum grid (typically, N ≈ 10 5 − 10 6 ).To obtain the matrix A, we time-step the rt-BTE in a sampling window consisting of M time steps (using the Perturbo code [38]) and then we form two matrices X 1 and X 2 relating the populations at times t and t + ∆t.The populations f k (t) from t 1 to t M −1 are stacked column-wise in the matrix X 1 , with column i corresponding to time t i and containing the populations f k (t i ) for all k-points and bands.The populations from t 2 to t M are similarly stacked columnwise in the second matrix X 2 .According to equation 2, these matrices are related by X 2 = AX 1 , but computing A naively from the pseudoinverse of X 1 has a prohobitive cost due to the large size N of the k-point grid.To circumvent this problem, in DMD one first performs a truncated singular value decomposition (SVD) [43,44] of the X 1 matrix: where Σ ∈ R N ×M is a matrix with diagonal entries equal to the singular values σ j arranged in decreasing order, while U ∈ C N ×N and V ∈ C M ×M are matrices collecting the mutually orthogonal singular vectors [45].(Above, V † indicates the Hermitian conjugate of V.) Because X 1 contains the time-dependent populations, this SVD procedure can single out the main patterns in the momentum-space dynamics.Here, we keep only the first r singular values (typically, r ≈ 10) to restrict the solution space to the leading r momentum-space modes, and then project the matrix A onto this reduced rdimensional space.This procedure provides the matrix Ã, with reduced size r × r, which can be diagonalized straightforwardly to obtain the dominant DMD modes.Using this procedure, the populations at future times t > t M are predicted − that is, obtained without explicit solution of the rt-BTE − using where ϕ l k are the momentum-space DMD modes obtained from the matrix Ã, and ω DMD l and b l are their frequencies and amplitudes (see Methods for detailed derivations).
We summarize the main steps of this DMD procedure, which are illustrated in Fig. 1: 1. Simulate the rt-BTE dynamics for the first M steps and construct the matrices X 1 and X 2 ; 2. Perform SVD on X 1 to find the matrix Ã in the reduced r-dimensional space; 3. Diagonalize Ã to find the DMD modes ϕ l k and their frequencies ω DMD l , with l = 1 . . .r; 4. Obtain the mode amplitudes b l from the initial condition f k (t 1 ); 5. Predict the dynamics for t > t M using equation 4.
A key parameter is the duration of the sampling window (t M ) required for accurate DMD extrapolation of the dynamics beyond t M .As the computational cost of DMD is negligible, the size of the sampling window, during which the rt-BTE is solved by explicit time-stepping, determines the computational cost of the entire workflow.

C. High-field electron dynamics
We employ our DMD-based approach to simulate timedomain electron dynamics in an applied electric field in the presence of e-ph collisions.We recently demonstrated similar calculations using the rt-BTE without the aid of data-driven techniques [31].Here we use this case study to explore the accuracy and efficiency of our rt-BTE plus DMD approach as well as find optimal values for the sampling window and analyze the momentumspace DMD modes.Our calculations focus on electrons in GaAs, where the conduction band has three sets of low-energy valleys, at Γ and near L and X in order of increasing energy [46] (see the inset in Fig. 2a).Upon applying an electric field, the electrons are accelerated to higher band energies while they also transfer part of that excess energy to the lattice via e-ph collisions.These competing mechanisms lead to a steady-state electronic distribution which is typically reached on a picosecond to nanosecond time scale.
Our simulations begin with electrons in thermal equilibrium with the lattice at 300 K. We apply a constant electric field E and time step the electron populations until they reach the steady state distribution, f E k , from which we compute the mean drift velocity, v(E), a quantity routinely measured in experiments [47][48][49].Repeating this procedure for multiple field values allows us to construct the full drift velocity versus electric field curve in a material, starting from linear response at low field to velocity saturation at high field [31].
Figure 2a shows the time-dependent populations in four regions of the Brillouin zone following the application of a high field (5 kVcm −1 ).Electrons initially occupying the Γ-valley scatter to the higher-energy Land X-valleys.As a result, the electron populations in the Γ-valley decrease, with a corresponding increase in L-and X-valley populations.In regions of momentum space between the Γ-and L-valleys the populations peak at intermediate times and then relax to lower values.This dynamics is nontrivial because the populations evolve differently in different momentum-space regions, making accurate predictions challenging.Our DMD approach can learn the dominant modes governing this intricate dynamics and extrapolate the time-dependent populations well beyond the sampling window.Remarkably, we find that a short sampling window − 400 fs to 2 ps out of a total simulation time of 12.5 ps − is sufficient to extrapolate the dynamics all the way to steady state, with rt-BTE and DMD trajectories in nearly exact agreement outside the sampling window (Fig. 2a).This accuracy extends to the entire set of ∼10 5 k-points considered in our simulations.
The ability to learn key temporal momentum-space patterns is a consequence of the relatively rapid decay of the singular values of the X 1 matrix used for learning the dynamics in the sampling window (Fig. 2b).This decay becomes slower as the sampling window increases, but it remains significant even for the longest sampling window of 2 ps used here (note the log scale in the plot).In turn, the singular value decay enables a striking dimensionality reduction, with DMD employing only r ≈ 10 modes to solve the dynamics as opposed to 10 5 populations f k and billions of e-ph scattering terms in the rt-BTE.
The choice of an ideal sampling window can rely on the appearance of specific DMD modes at steady state.Figure 2c shows the DMD mode frequencies ω DMD in the complex plane, where the imaginary part of ω DMD corresponds to the decay rate of a given mode and the real part gives its oscillation frequency.The populations f k (t) are real-valued and are written as a summation of complex exponentials in equation 4. Therefore, physically meaningful results are possible only when Re(ω DMD ) = 0 (modes 1, 2) or when ω DMD appear as complex conjugate pairs (modes 3 − 10).Describing the steady state is particularly important in our simulations.In DMD, all modes with a non-zero imaginary frequency vanish in the long time limit, with only one mode surviving at steady state (mode 1 in Fig. 2c).As the sampling window increases, the imaginary frequency of this mode goes to zero, providing the correct steady state behavior.This analysis allows us to find the minimal sampling window required for accurate steady-state results by monitoring the zero-frequency mode.
The DMD eigenvector of the zero-frequency mode (mode 1 in Fig. 2d) determines the steady-state electron distribution ϕ 1 k = f k (t → ∞), while the other modes control the transient dynamics.For example, mode 2 governs electron scattering from the Γ-to the Land X-valleys, and higher modes appearing as conjugate pairs exhibit oscillating trends in energy (modes 3−10 in Fig. 2d).Converging the zero-frequency mode allows us to compute the steady-state drift velocity more efficiently.Figure 2e shows that a sampling window of 1.7 ps (170 snapshots) provides a drift velocity nearly identical to the full rt-BTE calculation, which requires much longer simulation times of up to 12.5 ps (1250 snapshots).On this basis we conclude that DMD needs only ∼10% of the dynamics data for accurate steady-state predictions.

D. Velocity-field curves
We also employ DMD to accelerate calculations of entire velocity-field curves.This requires the drift velocity for a set of electric field values, and thus we adopt a modified workflow.Following our recent work [31], we gradually increase the electric field (black curve in Fig. 3a) and use the steady-state populations for a given field, f E k , as the initial condition for the next field value, E + ∆E, where the field increment ∆E is typically 100 − 200 Vcm −1 .As the applied field increases, the DMD frequencies and momentum-space modes change substantially.Therefore, for each new field value we repeat DMD learning in the initial stage of the simulation (see the DMD sampling regions shown as red rectangles in Fig. 3a).We then use DMD to predict the steady-state populations f E k and drift velocity for that field value, and use f E k as the initial condition for the next field value.
The velocity-field curves obtained with this approach are shown in Fig. 3b for GaAs and graphene and compared with rt-BTE results obtained without DMD.Using DMD lowers significantly the computational cost The DMD sampling window for each electric field is shown with a red rectangle, and the drift velocities outside this window are predicted with DMD.The steady-state drift velocities from DMD correspond to the plateaus for each field value, and agree with the reference drift velocities, obtained by explicitly time-stepping the rt-BTE until steady state, which are shown with white dots.b, Velocity-field curves in GaAs and graphene obtained from the rt-BTE (black and blue solid lines) and by combining the rt-BTE and DMD (red and orange dashed lines).
to obtain the full velocity-field curves, by a factor of 10.5 for GaAs and ∼16.5 for graphene, while fully preserving the accuracy.This efficiency is a consequence of DMD's ability to capture the dominant modes in the population dynamics using only a small number of snapshots, with a similar accuracy regardless of the electric field value.Our strategy of gradually increasing the electric field leads to an easier-to-extrapolate dynamics compared to the abrupt application of a strong field.With our approach, the DMD modes exhibit smooth trends in momentum space and the DMD frequencies are relatively small in magnitude, thus enabling accurate DMD predictions.

E. Excited electron relaxation
Next, we consider a different nonequilibrium dynamics where the material is initially prepared in an excited electronic state.This setting can be used, for example, to model the effect of an optical excitation with a laser pulse [13].Different from the high-field dynamics, in this case the long-time limit is known and we are primarily interested in the transient dynamics.Following the initial excitation, in the presence of e-ph interactions and without any external fields, the electrons relax to a thermal equilibrium Fermi-Dirac distribution [50], f FD k , typically on a sub-picosecond time scale.This ultrafast dynamics can be modeled by time-stepping the rt-BTE until reaching the equilibrium Fermi-Dirac distribution.Using this approach, our previous work has shown that electrons relax to the band edge significantly slower than holes in GaN semiconductor, with implications for optoelectronic devices [15].
Following that work, we model an excited state in GaN by placing the electrons ∼1 eV above the conduction band edge, and then obtain the time-dependent electron populations by solving the rt-BTE (see Fig. 4a,b).We employ DMD to predict this transient dynamics, and find large errors when using a short time window of up to ∼50 fs (solid red line in Fig. 4c).The correct steady state and transient dynamics are obtained by increasing the sampling window to 200 fs (dashed orange line in Fig. 4c).Our analysis of the DMD frequencies shows that the zero-frequency mode describing thermal equilibrium in the long-time limit appears when the sampling window reaches 100 fs (see the arrow in Fig. 4d) and fully converges for a ∼200 fs sampling window.The need for such a long sampling window relative to the total duration of the dynamics (400 fs) makes DMD ineffective.
To address this issue and more efficiently study transient dynamics with DMD, we formulate a different learning procedure that incorporates knowledge of the equilibrium state.We focus on the difference between the transient and equilibrium populations, δf k (t) = f k (t) − f FD k , as opposed to just f k (t) as we did in the high-field example.After predicting δf k (t) with DMD, we obtain the time-dependent populations f k (t) by adding back the f FD k term.As δf k vanishes in the long-time limit (Fig. 4b), the zero-frequency DMD mode is missing when computing δf k (Fig. 4e); all other DMD frequencies associated with δf k are similar to those for f k (t) (Fig. 4d-e).We find that the DMD method based on δf k is far more effective and requires a significantly shorter sampling window for accurate DMD predictions − using a 50 fs sampling window, we achieve results similar to DMD for f k (t) with a four times longer (200 fs) window (Fig. 4c).
With this improved DMD approach, using a sampling window of only ∼12% of the total simulation time allows us to accurately predict the average electron relaxation rate in GaN, with a DMD computed value of 5.23 eVfs −1 in close agreement (within 0.8%) with the rt-BTE result.This result demonstrates that our DMD approach can predict excited electron relaxation with a high accuracy.

III. DISCUSSION
The DMD approach introduced here is highly efficient − different from the rt-BTE, which requires a supercomputer with extensive memory and CPU resources, our DMD workflow can be performed straightforwardly on a laptop.The most demanding step is carrying out truncated SVD on the X 1 matrix, but for comparison this step requires lower computational resources than even just a single rt-BTE time step.This remarkable speed-up is achieved by reducing the dimensionality of the rt-BTE dynamics and is linked to the shape of the X 1 matrix.The rt-BTE employs a large number of k-points (about 10 5 − 10 6 ), which equals the number of rows of the matrix X 1 , and a significantly smaller number of snapshots in the DMD sampling window, typically ∼100 time steps, which sets the number of columns in X 1 .Following truncated SVD, the size of the problem is reduced to (at most) the number of snapshots and is typically of order 50−100, and thus smaller by orders of magnitude compared to the original rt-BTE.(Note that one could use the entire set of singular values, but here we prefer using only ∼10 singular values to prevent numerical instabilities [45]).
This efficiency allows us to evaluate the accuracy of DMD on the fly, halting explicit time-stepping of the rt-BTE when the DMD steady state or transient dynamics are fully converged.In addition, our approach addresses the key challenge of storing the rt-BTE populations, which are needed only in the sampling window in DMD, as opposed to the full dynamics.This is a critical improvement because in conventional rt-BTE simulations one needs to store the populations f k (t) on dense momentum grids for thousands of time steps, resulting in terabytes of data.In contrast, after carrying out SVD in the sampling window, DMD stores only a handful of complex frequencies and momentum-space modes, using which the dynamics can be reconstructed for the entire simulation.

IV. CONCLUSION
In summary, we have introduced a data-driven approach based on DMD to accelerate first-principles calculations of nonequilibrium electron dynamics in materials.Our method speeds-up the solution of the time-dependent Boltzmann equation with electron collisions computed from first principles.We have shown that DMD can capture dominant modes governing the microscopic dynamics, enabling accurate predictions of the steady-state properties such as the drift velocity as well as transient processes such as electron relaxation and equilibration.In both steady-state and transient nonequilibrium calculations, DMD requires explicit timestepping of the rt-BTE in a time window of only ∼10% of the full simulation, after which the dynamics is extrapolated from the DMD modes with negligible computational cost.This DMD workflow preserves the accuracy while requiring far more modest computational resources than full rt-BTE simulations.
These advances are broadly relevant to studying nonequilibrium quantum dynamics of elementary excitations.For example, in future work our data-driven approach will be adapted to study phonon dynamics, which requires costly computations of phonon-phonon scattering, as well as exciton dynamics using a bosonic rt-BTE formalism.
V. METHODS

A. Computing DMD modes and frequencies
Let us describe in more detail the calculation of DMD modes and frequencies.We start from the snapshots f k (t) evaluated explicitly with the rt-BTE in the sampling window t 1 < t < t M , and then apply the SVD procedure to the matrix X 1 (see equation 3).As shown in Fig. 2b, we find that the singular values σ j decay rapidly.Keeping only the largest r ≈ 10 singular values, we write the SVD of X 1 as where we defined the economy-sized matrices in the rdimensional subspace [45] as Σ = Σ(1 : r, 1 : r), Ũ = U(1 : N, 1 : r), Ṽ = V(1 : M, 1 : r).This way, the approximate pseudo-inverse of the matrix X 1 , denoted as X 1 + , can be obtained with little effort as Ṽ Σ−1 Ũ † .Then the matrix A relating the snapshot matrices via X 2 = AX 1 can be written as Due to its large N × N size (here, N ≈ 10 5 is the number of k-points), diagonalizing A is computationally expensive.In DMD, a key step is rewriting this matrix in the reduced r-dimensional space: allowing for straightforward eigenvalue decomposition: ÃW = WΛ, where the matrix W contains the eigenvectors of Ã and the eigenvalues Λ = diag{λ l } are common to both matrices Ã and A [51].The DMD modes, stacked column-wise in the matrix The DMD frequency of mode l is obtained from the corresponding eigenvalue λ l using equation 8, and thus the mode amplitude vector b is obtained from the pseudo-inverse of the DMD mode matrix Φ: This approach provides the DMD modes ϕ l k , frequencies ω DMD l , and mode amplitudes b l , and thus all the quantities needed for DMD prediction of the dynamics outside the sampling window (t > t M ) using equation 4.

B. Electron-phonon scattering from first principles
Our first-principles calculations of e-ph scattering employ an established workflow, which is summarized here and described in more detail in Ref. [38].The electronic wave functions and band energies are obtained from plane-wave DFT calculations with the Quantum Espresso code [52] using the local density approximation [53] and norm-conserving pseudopotentials [54].The electronic quasiparticle band structure is refined using GW calculations carried out with the YAMBO code [55].This step improves the agreement with experiment of the electron effective masses and relative valley energies, which are essential for precise calculations of high-field dynamics [31] and excited electron relaxation [15].
The phonon dispersion and e-ph perturbation potentials are obtained from density-functional perturbation theory (DFPT), where lattice vibrations and their coupling with electrons are treated as perturbations to the ground-state electron density [52].The e-ph interactions are described by the matrix elements g mnν (k, q) = ℏ 2ω νq 1/2 ψ mk+q ∆ νq V KS ψ nk , (13) which are the probability amplitudes to scatter from an initial electronic state |ψ nk ⟩, with band n and momentum k, to a final state |ψ mk+q ⟩ by absorbing or emitting a phonon with mode ν, momentum q, and frequency ω νq .The term ∆ νq V KS is the e-ph perturbation potential induced by the phonon mode and is defined in Ref. [38].
To study nonequilibrium dynamics with the rt-BTE, the electrons, phonons, and e-ph scattering are described on dense kand q-point momentum grids.Obtaining the e-ph matrix elements on such grids directly from DFPT is computationally prohibitive.Therefore, we first compute g mnν (k, q) on coarse kand q-point grids, and then interpolate these quantities to significantly finer grids using Wannier-Fourier interpolation with Wannier functions generated from Wannier90 [56].Finally, the e-ph scattering integral employed in equation 1 is defined as where N q is the number of q-points, ε nk and ε mk+q are the band energies of the initial and final electronic states, and the Dirac delta functions expressing energy conservation are implemented as Gaussians with a small (∼5 meV) broadening.Above, F abs and F em are phonon absorption and emission terms, whose explicit expressions are given in Ref. [38].The Wannier interpolation, scattering integral computation, and the rt-BTE ultrafast dynamics are implemented in our Perturbo opensource package [38].

FIG. 1 .
FIG.1.Workflow of DMD plus rt-BTE calculations.The first M steps of the dynamics, which make up the sampling window for DMD learning, are simulated by solving the rt-BTE.The resulting populations f k (t) are stacked in the X1 and X2 matrices with a relative shift of one time step.The dynamics at later times t > tM is predicted with DMD using the r leading modes obtained by SVD of the matrix X1 and diagonalization of the matrix Ã = Ũ † A Ũ.

FIG. 2 .
FIG. 2. DMD simulations of electrons in GaAs in an applied electric field.a, Time-dependent electron populations for four electron momenta.The gray region indicates the shortest sampling window (0.4 ps) and the red vertical line the longest sampling window we tested (2 ps).Solid black lines show the rt-BTE results and orange dashed lines the DMD predictions obtained using the longest sampling window.The inset is a schematic of the low-energy band structure of GaAs showing the Γ-and higher energy L-and X-valleys.b, Singular values of the X1 matrix, with the ten largest singular values used in our DMD calculations separated by a vertical line.c, DMD frequencies plotted in the complex plane and shown as circles with radii proportional to the DMD mode amplitudes b l .In panels a, b, and c, the colors indicate the duration of the DMD sampling window according to the legend given in (c).d, DMD momentum-space modes ϕ l k , multiplied by the corresponding amplitudes b l , given as a function of energy.The initial state f k (t1) is shown with a dashed line.e, Convergence of the steady-state drift velocity with respect to duration of the DMD sampling window.The rt-BTE value is shown for reference as a dashed line.

FIG. 3 .
FIG.3.Velocity-field curves from DMD. a, Transient drift velocity in GaAs computed as a function of time (black curve).The external electric field is increased step-wise in the simulation (see the field values given above the plot).The DMD sampling window for each electric field is shown with a red rectangle, and the drift velocities outside this window are predicted with DMD.The steady-state drift velocities from DMD correspond to the plateaus for each field value, and agree with the reference drift velocities, obtained by explicitly time-stepping the rt-BTE until steady state, which are shown with white dots.b, Velocity-field curves in GaAs and graphene obtained from the rt-BTE (black and blue solid lines) and by combining the rt-BTE and DMD (red and orange dashed lines).

FIG. 4 .
FIG.4.dynamics in GaN using DMD.a, Energy dependence of the momentum-averaged electron populations, and b, difference between the time-dependent and equilibrium populations in GaN.In both panels, the simulation time is color-coded using sepia for the initial excited state and purple for the equilibrium state.The energy zero is set to the conduction band minimum.c, DMD error on the electron populations, computed as the root-meansquare difference between the reference values from rt-BTE and those obtained from DMD. Results are shown for different sampling windows, given in the legend, and for the two schemes where DMD is applied to f k (t) or alternatively to δf k (t) = f k (t) − f FD k .d, e, DMD frequencies on the complex plane, respectively for f k (t) and f k (t) − f FD k , with the duration of the sampling window color-coded.
where ∆t is the simulation time step.The mode amplitudes b = b 1 b 2 • • • b r ∈ C r are obtained from the initial condition.Setting t = 0 in equation 4, we get