# Theoretical model of passive mode-locking in terahertz quantum cascade lasers with distributed saturable absorbers

Lukas Seitner, Johannes Popp, Michael Haider, Sukhdeep S Dhillon, Miriam S Vitiello, Christian Jirauschek

► To cite this version:

Lukas Seitner, Johannes Popp, Michael Haider, Sukhdeep S Dhillon, Miriam S Vitiello, et al.. Theoretical model of passive mode-locking in terahertz quantum cascade lasers with distributed saturable absorbers. Nanophotonics, 2024, 13 (10), pp.1823 - 1834. 10.1515/nanoph-2023-0657. hal-04732879

HAL Id: hal-04732879

https://hal.science/hal-04732879v1

Submitted on 11 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Research Article

Lukas Seitner*, Johannes Popp, Michael Haider, Sukhdeep S. Dhillon, Miriam S. Vitiello and Christian Jirauschek*

# Theoretical model of passive mode-locking in terahertz quantum cascade lasers with distributed saturable absorbers

https://doi.org/10.1515/nanoph-2023-0657

Received October 4, 2023; accepted January 19, 2024;

published online February 1, 2024

**Abstract:** In research and engineering, short laser pulses are fundamental for metrology and communication. The generation of pulses by passive mode-locking is especially desirable due to the compact setup dimensions, without the need for active modulation requiring dedicated external circuitry. However, well-established models do not cover regular self-pulsing in gain media that recover faster than the cavity round trip time. For quantum cascade lasers (QCLs), this marked a significant limitation in their operation, as they exhibit picosecond gain dynamics associated with intersubband transitions. We present a model that gives detailed insights into the pulse dynamics of the first passively mode-locked QCL that was recently demonstrated. The presence of an incoherent saturable absorber, exemplarily realized by multilayer graphene distributed along the cavity, drives the laser into a pulsed state by exhibiting a similarly fast recovery time as the gain medium. This previously unstudied state of laser operation reveals a remarkable response of the gain medium on unevenly distributed intracavity intensity. We show that in presence of strong

spatial hole burning in the laser gain medium, the pulse stabilizes itself by suppressing counter-propagating light and getting shortened again at the cavity facets. Finally, we study the robustness of passive mode-locking with respect to the saturable absorber properties and identify strategies for generating even shorter pulses. The obtained results may also have implications for other nanostructured mode-locked laser sources, for example, based on quantum dots.

**Keywords:** terahertz; quantum cascade laser; graphene; frequency comb; passive mode-locking; Maxwell-Bloch

## 1 Introduction

The operation of lasers in a mode-locked state is highly desirable in many applications, such as metrology and sensing [1], [2]. In the spectral domain, mode-locking is associated with the formation of comb spectra consisting of equidistant discrete lines. The realization of laser frequency combs was awarded the Nobel Prize in physics, indicating its importance in present and future applications [2], [3]. Fundamental research of light-matter interactions can be performed for many different timescales, down to femtoseconds, using ultrashort laser pulses [4]. Specifically, in the terahertz (THz) and mid-infrared (mid-IR) portion of the electromagnetic spectrum, the dynamics of gas molecules may be explored as they exhibit roto-vibrational modes at these frequencies. The corresponding optical transitions can be used for molecular fingerprint detection in sensing and spectroscopy [5]. While spectroscopic platforms have been realized in this fingerprint region (e.g., [6], [7]), the absence of compact and inexpensive short pulse emitters hinders the implementation of short pulse applications [8]. Recently, a passively mode-locked THz quantum cascade laser (QCL) was demonstrated in a compact configuration, which may fill this gap [9]. Here, we present a model that gives a detailed understanding of the mode-locking dynamics in this device, as is required for an improved understanding

<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>

<footnote>*Corresponding authors: Lukas Seitner, TUM School of Computation, Information and Technology, Technical University of Munich (TUM), D-85748 Garching, Germany, E-mail: lukas.seitner@tum.de.</footnote>

<footnote>https://orcid.org/0000-0002-0985-8594; and Christian Jirauschek, TUM School of Computation, Information and Technology, Technical University of Munich (TUM), D-85748 Garching, Germany; and TUM Center for Quantum Engineering (ZQE), Technical University of Munich (TUM), D-85748 Garching, Germany, E-mail: jirauschek@tum.de. https://orcid.org/0000-0003-0785-5530</footnote>

<footnote>Johannes Popp and Michael Haider, TUM School of Computation, Information and Technology, Technical University of Munich (TUM), D-85748 Garching, Germany. https://orcid.org/0000-0003-1745-4888 (J. Popp). https://orcid.org/0000-0002-5164-432X (M. Haider)</footnote>

<footnote>Sukhdeep S. Dhillon, Laboratoire de Physique de l'Ecole Normale Supérieure, ENS, Université PSL, CNRS, Sorbonne Université, Université de Paris, Paris, France.</footnote>

<footnote>Miriam S. Vitiello, NEST, CNR - Istituto Nanoscienze and Scuola Normale Superiore, Piazza San Silvestro 12, 56127, Pisa, Italy.</footnote>
of the underlying mechanisms and a targeted design and optimization of such devices. The obtained results may also have implications for the design of other passively mode-locked semiconductor lasers with short gain recovery times, for example, based on quantum dots (QDs).

In general, amplitude-modulated mode-locking can be achieved by active injection of an optical or electrical signal, or passively by a suitable modification of the laser cavity [1]. Passive mode-locking does not require dedicated circuitry and allows for more compact setups, thus having advantages for practical applications. Furthermore, the minimum achievable pulse duration is inherently shorter than that achievable via active mode-locking [1]. In general, passive mode-locking can be realized through saturable absorbers [10], [11] or nonlinear phase shifts [12], [13]. Today's standard theory of passive mode-locking explains the dynamics of pulse formation very well when the gain recovers more slowly than the cavity round trip time [1]. However, the model does not cover self-starting passive mode-locking when the gain dynamics become too fast. QCLs extend the spectral regime covered by interband semiconductor lasers to mid-infrared and THz frequencies and exhibit gain recovery times significantly faster than the cavity round trip time [14], [15]. Thus, passive mode-locking of QCLs was assumed to be rather difficult, and research focused mainly on frequency-modulated (FM) combs in Fabry-Pérot cavities and soliton generation in ring cavities, e.g., [16]–[20]. The possibility of comb generation in the low THz region by difference frequency generation of mid-infrared radiation inside a single QCL has been explored in experiment and simulation [21]–[23]. But also direct generation of THz frequency combs can be achieved by dispersion engineering, even though the resulting waveform typically features both frequency and amplitude modulation rather than single isolated pulses [24]. Several theoretical works have been dealing with the possibility of amplitude-modulated mode-locking in fast gain media, based on self-induced transparency [25], a combination of active and absorbing regions [26], [27], coherent and memory effects [27]–[29], or optical injection [30]. However, the generation of passively mode-locked pulses has been experimentally achieved only very recently [9]. In this device, multilayer graphene (MLG) stripes embedded alongside the top contact serve as a distributed saturable absorber for the THz radiation and drive the laser in the pulsed state. It has been proven experimentally that THz saturable absorption in graphene depends on the number of layers ($n_{\text{layer}}$) and that transparency modulation increases with $n_{\text{layer}}$ [31], [32]. Hence, MLG is employed to enhance the saturable absorption effect, which is known to be the key mechanism for passive mode-locking. Apart

from that, there is also a technical limitation in the use of single-layer graphene: its clean adhesion onto the top surface stripes of a processed QCL is worse than in the case of MLG, which conversely ensures good adhesion, flatness, and preserves a clean process. In a previous report [9], we presented numerical results based on a Maxwell-density matrix approach, considering nine quantized levels per period in the active region. This model showed very good agreement with experimental data but, due to its complexity, did not provide an intuitive understanding of the mode-locking dynamics [9]. Starting from this framework, we develop a simplified description, which still contains the relevant dynamics but offers an intuitive understanding of the pulse formation and stabilization. This allows us to pin down the necessary properties of absorber and gain for pulse formation. Since the Maxwell-Bloch equations are widely used for describing nano-optoelectronic devices [33], our model provides a basis for investigating related effects in other types of semiconductor lasers.

## 2 Theoretical model

One major conclusion drawn in the seminal paper by Haus is that the “model of a fast saturable absorber and a fast laser medium does not give self-starting mode-locking solutions” [34]. Notably, for intermediate situations with gain recovery times on a comparable timescale as the round trip time, self-starting mode-locking is not excluded by the Haus theory, but this case is less amenable to analytical solutions and has less been studied in literature (see, e.g., [1], [28], [29]). However, for gain recovery times much shorter than the round trip time the conditions under which self-starting passive mode-locking can still be obtained become increasingly restricted (see, e.g., [26]). In the limiting case of fully recovered gain and loss after a round trip, the existence of a net gain window exclusively at the position of the pulse is not compatible with a self-starting scenario, requiring that the unsaturated gain exceeds the unsaturated loss. Thus, this section focuses on additional effects due to counter-propagation and associated spatial hole burning, not included in the Haus approach but relevant especially in semiconductor lasers, which can lead to additional pulse stabilization.

First, we introduce the optical propagation equation, which is coupled to the absorber model. Then, we describe the active region dynamics based on the Bloch equations with spatial hole burning (SHB) and extract an analytical formula for the direction-dependent gain.

The overall model implemented in the present work is summarized in Figure 1.
<figure><img src="image_5.png" /><figcaption>**Figure 1:** Visualization of the semi-classical, dynamical model, which consists of three domains and their interactions. The electron dynamics of a fast gain medium are described by the Bloch equations. A possible realization is the quantum-well heterostructure of a QCL active region. The reduced two-level model accounts for coherent optical transitions and incoherent scattering between adjacent periods. The polarization and electric field terms introduce coupling between the active region description and the optical propagation equation, which is derived from Maxwell's equations. A saturable absorber model accounting for the distributed multilayer graphene completes the system description. It interacts with the optical field via an intensity-dependent loss.</figcaption></figure>

## 2.1 Optical propagation in presence of saturable absorption

The time evolution of the electric field is governed by a propagation equation, which is derived from Maxwell's equations as [1], [35], [36]

$$v_g^{-1} \partial_t E^\pm = \mp \partial_x E^\pm + f^\pm - lE^\pm + \left(D_f - i\frac{\beta_2}{2}\right) \partial_t^2 E^\pm. \quad (1)$$

The electric field is in the slowly varying envelope approximation (SVEA) described as the sum of two counter-propagating components $E(x, t) = \sum_{\pm} \frac{1}{2} [E^\pm \exp(\pm i\beta_0 x - i\omega_c t) + \text{c.c.}]$, where c.c. denotes the complex conjugate, $\beta_0$ is the propagation constant, $\omega_c$ represents the center frequency, and $v_g$ denotes the group velocity. The term $D_f$ introduces the effect of gain bandwidth, while $\beta_2$ models the background material's group velocity dispersion. The field components of Eq. (1) interact with the saturable absorber via the loss term $lE^\pm$, and with the gain medium via the intersubband polarization $f^\pm$. These two driving forces of the system need to be balanced in order to achieve the pulsed state. Notably, this propagation equation holds significant similarity to the master equation derived by Haus [1], regarding the treatment of saturable loss, dispersion, and gain bandwidth $D_f$. However, it does not explicitly contain a nonlinearity term [1], as, e.g., considered for the modeling

of QCL FM combs [36]–[38], such that all nonlinear effects in the here-discussed system arise from the standard Bloch equations with SHB [39] and the addition of the saturable absorber.

Due to the outstandingly high electron mobility and ballistic transport properties in cryogenic graphene [40], [41], we assume absence of SHB in the saturable absorber. Furthermore, because of the peculiar band structure [42] and complex optoelectronic dynamics [43], [44], we choose a compact phenomenological saturable absorber model for the graphene at THz frequencies [9], [44]–[46]. For finite recovery times, saturable absorption is typically described by [35]:

$$\partial_t a_s = - \frac{I}{\tau_a I_{s,a}} a_s - \frac{a_s - a_1}{\tau_a}. \quad (2)$$

The saturable loss then enters the propagation equation via the relation

$$l(x, t) = \frac{1}{2} [a_0 + a_s(x, t)]. \quad (3)$$

As can be seen, the overall loss is split into a constant waveguide and material contribution $a_0$, and a saturable part $a_s$. Due to the distribution of graphene over the whole cavity, this part is a function of space and time $a_s(x, t)$ and has a maximum value of $a_1$. Therefore, Eq. (2) accounts for loss saturation with increasing intensity $I(x, t)$ and recovery with a characteristic time $\tau_a$. The remaining parameter in this model is the effective saturation intensity $I_{s,a}$, defined
as the value of the intracavity intensity that saturates the graphene loss to $a_s = a_1/2$.

## 2.2 Modeling of electron dynamics

To describe the QCL gain medium, we aim for a model that is simple enough to give intuition but still accounts for the physical mechanisms that enable passive mode-locking. From a very detailed multilevel density matrix approach [33], [47], as used in [9], we construct a more compact two-level Bloch model [33], [39], [48], [49]. The corresponding parameters are determined by preserving the crucial properties of the resulting gain, i.e., maximum, recovery time, and bandwidth. Thus, no parameter fitting is needed. The state of the two-level system is determined by the coherence $\eta$, inversion $w$, and inversion grating $w^+$ [39]. For bidirectional optical propagation, as present in Fabry-Pérot cavities, this inversion grating needs to be considered, as it accounts for optical interference effects like SHB. This effect is known to play a major role in multimode and comb operation of QCLs and other nanostructured lasers and is, therefore, explicitly considered in our model [16], [36], [38], [39], [50], [51]. In homogeneously broadened laser systems, the decay of quantum coherence between the laser levels is typically the fastest process in the system. In a THz QCL, it is usually about 10 times faster than the gain recovery time. By assuming it to be infinitely fast, i.e., performing the adiabatic elimination of $\eta$ [52], the Bloch equation system can be simplified, allowing for a more intuitive understanding. After restructuring the equations and disregarding minor contributions, they read as:

$$ \eta_{21}^{\pm} = -\frac{i}{2\hbar\gamma_2} d_{21} [wE^{\pm} + w^{\pm}E^{\mp}], \quad (4) $$

$$ \partial_t w = -\frac{d_{21}^2}{\hbar^2 \gamma_2} w (|E^+|^2 + |E^-|^2) - \gamma_1 (w - w_{eq}), \quad (5) $$

$$ \partial_t w^+ = -\frac{d_{21}^2}{\hbar^2 \gamma_2} w E^+(E^-)^* - \gamma_1' w^+. \quad (6) $$

Here $w^- = (w^+)^*$, $d_{21}$ is the optical dipole moment, $\gamma_1$ denotes the gain recovery rate, $\gamma_2$ is the dephasing rate, and $w_{eq}$ represents the equilibrium inversion. Further, $\gamma_1'$ is the recovery time of the inversion grating, which is affected by charge carrier diffusion,

$$ \gamma_1' = \gamma_1 + 4\rho_0^2 D, \quad (7) $$

were $D$ denotes the diffusion constant. A more detailed derivation of these equations is given in the Supplementary Information Section I.B. The resulting polarization of this system, which enters the electric field propagation equation, is then given by

$$ f^{\pm} = i\alpha\eta_{21}^{\pm}, \quad (8) $$

with

$$ \alpha = \frac{\Gamma n_{3D} \omega_c d_{21}}{c_0 \varepsilon_0 n_0}. \quad (9) $$

Here, $\Gamma$ denotes the overlap factor, $n_{3D}$ is the space averaged doping density, and $c_0$, $\varepsilon_0$, and $n_0$ are the vacuum speed of light, permittivity, and material refractive index, respectively. We notice from Eqs. (4)–(9) that in presence of the inversion grating $w^+$, the gain is not only saturated by an intensity $\propto (|E^+|^2 + |E^-|^2)$ but also an interference term $E^+(E^-)^*$ is present. In fact, when setting $w^+ = 0$ in, then Eq. (4), Eqs. (4) and (5) correspond to an implementation of the gain into the Haus model in analogy to Eq. (2) (compare to Supplementary Information Section I.D). Therefore, our model with SHB can be seen as a direct extension to the Haus theory.

As a consequence of the adiabatic elimination, gain curvature and dispersion get neglected. However, since a finite gain bandwidth is necessary for stable mode-locked solutions, it is reintroduced via $D_f$ in the propagation equation (1) [1], [35]. The derivation of $D_f$ is outlined in the following subsection, along with further analytical considerations, that are crucial to describe and understand the pulse dynamics. Furthermore, by adjusting $v_g$ and $\beta_2$ in Eq. (1), gain dispersion may be considered, but this effect turns out to be of minor importance in our case.

## 2.3 Direction-dependent gain model

Our approach to calculate the gain is based on taking the field propagation equation, Eq. (1), and equating the inter-subband polarization to an optical amplification, which yields the forward and backward gain coefficient

$$ g^{\pm} = \pm 2\Re \left( \frac{f^{\pm}}{E^{\pm}} \right). \quad (10) $$

This direction-dependent gain means that left and right traveling fields might experience different amplification strengths. This approach may be counterintuitive at first sight, but in the field of nonlinear optics, counter-propagating effects can significantly influence the system dynamics [53]. For the implementation of the gain bandwidth in Eq. (1), we use the relation

$$ D_f = \frac{g_s}{2\gamma_2^2} \quad (11) $$

with a direction-averaged gain value

$$ g_s = \frac{g^+ I^+ + g^- I^-}{I^+ + I^-}, \quad (12) $$

where $I^\pm = \frac{1}{2}\varepsilon_0 c_0 n_0 |E^\pm|^2$ are the intensities of each field component.
Recent theoretical works have investigated the role of such counter-propagating interactions for different operating regimes of QCLs [20], [38], [36]. In the case of monochromatic excitation ($\partial_t w = 0$), we can derive an analytic expression for the direction-dependent gain from the adiabatic Bloch equations, Eqs. (4)–(6) (see Supplementary Information Section I.C):

$$g^{\pm} = \alpha \frac{d_{21}}{\hbar \gamma_2} \left( 1 - \frac{2I^{\mp}}{2I_{s,g} \gamma'_1 / \gamma_1 + I^+ + I^-} \right) w. \quad (13)$$

We observe that for two fields $E^{\pm}$ with different amplitudes, the stronger field component experiences a higher gain than the weaker counter-propagating field. This shows that the field with higher amplitude couples more efficiently with the gain medium. In a standard linear cavity, this effect is usually not essential, as facet reflections and irregular intensity fluctuations due to SHB dominate the dynamics. By adding the distributed saturable absorber, uneven intensity distributions arise, and the effect becomes more pronounced. For comparison, in a ring cavity, one field component completely decays in the absence of balancing mechanisms [19], [20], which is a direct result of the directive gain.

To get a better understanding of this phenomenon, we study two steady-state scenarios analytically, starting from expression (13). At first, we consider the case where the intensities in both field directions are equal ($I^{-} = I^{+}$), which happens at the facet of the cavity if we assume full reflectivity. This condition yields the maximum optical interference between the counter-propagating fields and effectively reduces the saturation intensity by the factor (see Supplementary Information Section I.C)

$$I_{s,g}|_{I^{-}=I^{+}} = \frac{2\gamma'_{1}}{2\gamma'_{1} + \gamma_{1}} I_{s,g}. \quad (14)$$

Assuming small diffusion and thus close-to-maximum SHB in the gain medium, we obtain $\gamma'_1 \approx \gamma_1$, and consequently $I_{s,g}$ is effectively reduced to 2/3 of its original value. For a suitable choice of parameters, the value of the gain can drop below the losses, even though exceeding it in an unsaturated state. In the case of a short pulse hitting the facet, the finite reaction times of gain and loss take effect, and a net loss window opens up at the trailing edge, leading to pulse shortening and suppression of subsequent fluctuations.

In the second case, we assume that one intensity direction is significantly stronger than the other one, $I^{\mp} \ll I^{\pm}$. This resembles the case where a high-intensity pulse travels along the cavity, and the counter-propagating intensity mainly consists of noise. Obviously, in this situation, no significant optical interference is present, and thus $I_{s,g}$ remains

at its initial value. However, since in this case (see Supplementary Information Section I.C)

$$g^{\mp}|_{I^{\mp} \ll I^{\pm}} = \left( 1 - \frac{2I^{\pm}}{2I_{s,g} \gamma'_{1} / \gamma_{1} + I^{\pm}} \right) g^{\pm}, \quad (15)$$

the gain gets reduced for the counter-propagating weak background. For example, if the pulse intensity is as strong as the saturation intensity, $I^{\pm} = I_{s,g}$, and we again assume negligible diffusion ($\gamma'_1 \approx \gamma_1$), then the resulting gain for the counter-traveling field would drop to $g^{\mp} = g^{\pm}/3$. In steady-state, the overall gain is clamped to the loss, meaning that the fainter field experiences a net loss.

For THz QCLs, $\gamma'_1 \approx \gamma_1$ is typically a valid approximation, since the comparably low frequency results in a small $\beta_0$ in Eq. (7), and the charge carrier diffusion is negligible at cryogenic temperatures. In semiconductor lasers that operate at higher frequencies and temperatures, like mid-IR QCLs, interband cascade (IC), or QD lasers, the relative impact of diffusion is increased. Nevertheless, SHB is still known to influence the dynamics of these lasers [39], [50], [51], [54], [55]. Thus, both previously described concepts are valid for these structures, and the impact of saturation intensity reduction and directional noise suppression can be quantitatively evaluated by inserting suitable parameters into Eqs. (14) and (15). Furthermore, the presented framework can be applied to saturable absorbers where a certain amount of SHB is present, to evaluate if the here-described pulse stabilizing mechanisms are still effective. This can be achieved by applying Eqs. (4)–(6), or respectively Eq. (13), also to the absorber, and evaluating the direction-dependent net gain $g^{\pm} - (a_0 + a_s^{\pm})$. For modeling passive mode-locking in QD lasers [56], the inclusion of SHB would require adapted Maxwell-Bloch equations [51], since the decay of the population grating in QD media is not well captured by the conventional diffusion model.

In the following, to reveal the full dynamics of the system, the coupled equations of gain, absorber, and intensity are studied numerically.

## 3 Numerical results and discussion

For our dynamic simulations of a mode-locked THz QCL, realistic parameter values extracted from experiment [9], literature [32], [44], [45], [57], and Monte-Carlo simulations [58] are used in the model. An overview of these values is given in Tables 1 and 2. Here, we describe the gain in terms of macroscopic parameters (maximum gain $g_0$, saturation intensity $I_{s,g}$, recovery time $1/\gamma_1$, coherence time $1/\gamma_2$). The microscopic contributions leading to these values (optical dipole moment, scattering rates, etc.) can be found in the
<table><caption>**Table 1:** Parameters of gain.</caption><tr><td>1/γ₁</td><td>1/γ₂</td><td>g₀</td><td>I<sub>s,g</sub></td></tr><tr><td>2.63 ps</td><td>0.167 ps</td><td>22.5 cm<sup>-1</sup></td><td>13 × 10<sup>7</sup> W m<sup>-2</sup></td></tr></table>

<table><caption>**Table 2:** Parameters of absorber.</caption><tr><td>τ<sub>a</sub></td><td>a<sub>0</sub></td><td>a<sub>1</sub></td><td>I<sub>s,a</sub></td></tr><tr><td>2.5 ps</td><td>11.8 cm<sup>-1</sup></td><td>10.5 cm<sup>-1</sup></td><td>5 × 10<sup>7</sup> W m<sup>-2</sup></td></tr></table>

Supplementary Information Section II. The losses ($a_0, a_1$) have been evaluated by dedicated electromagnetic simulations [9], and experiments show that the absorber recovers within $τ_a = 2$ ps to 3 ps [9], [32], [57], so we choose an intermediate value. Although the nonlinear optical properties of graphene have been successfully exploited to affect THz QCL operation [9], [59], [60], assessing its saturation intensity at THz frequencies appears to be complex, and different values can be found [32], [44]–[46]. Furthermore, the polarization direction of the electric field compared to the graphene surface has an influence on the saturation dynamics [44]. However, a polarization rotation possibly happening at the metal/graphene interface is not included in our one-dimensional model. Thus, we treat $I_{s,a}$ as the only fitting parameter in our model. Notably, the recovery times of gain and absorber are almost equal in the case considered here, and the gain maximum is slightly above the sum of background and saturable loss. A discussion on the robustness of the device against parameter variations will be given below and in the Supplementary Information.

## 3.1 Simulation of experimental setup

At first, we directly compare simulation results with time-resolved measurements, as given in [9]. In addition to the above stated parameters, we use realistic values for the mirror reflectivity ($R = 0.7$ [61]) and the diffusion constant in the gain medium ($D = 20 \text{ cm}^2 \text{ s}^{-1}$). The resulting time trace of the outcoupled intensity is given in Figure 2(a) and the corresponding spectra of the fields in Figure 2(b). A clear, background-free pulse operation with roughly 6 ps full-width-half-max (FWHM) duration is visible in the experiment and simulation, which is linked to a frequency comb in the spectral domain. This very good agreement supports the validity of the adiabatic simulation approach and also the reduction of the original nine-level system, which is self-consistently extracted from charge-carrier transport simulations [9], to an equivalent two-level model. The measured and simulated frequency combs, given in Figure 2(b),

<figure><img src="image_6.png" /><figcaption>**Figure 2:** Comparison of simulation results to experiment. (a) Time trace of the normalized intensity, outcoupled from the mode-locked laser. (b) Power spectrum of the electric field, normalized to the round-trip frequency with an offset of roughly 3 THz.</figcaption></figure>

exhibit around 10 lines in a 10 dB power range, with the experimental one having slightly more contributions at the higher frequency side. This small deviation may be due to the effect of after-pulsing in the measured time trace of Figure 2(a), which is not prominent in the simulations. When the integrated graphene absorber is removed (compare Supplementary Information Section III.A), the laser output becomes completely unlocked in both experiment and simulation. This straightforwardly leads to the conclusion that the graphene saturable absorber drives the QCL into the pulsed state.

In order to get the best visualization of the dynamics that keep a generated pulse stable (compare Section 2.3), we set the facet reflectivity to unity, while converting the outcoupling to a distributed loss [49], and neglect charge carrier diffusion in the quantum-well heterostructure. The resulting steady-state behavior of the intracavity intensity, gain, and loss is shown in Figure 3. In Figure 3(a), a pulse is located in the middle of the cavity, filling a significant part and traveling toward the right facet. The presence of this intensity peak leads to a dip in the distributed gain and loss. As both have fast recovery dynamics with nearly equal time constants, they have returned to their initial values at the cavity facets. The right-traveling intensity ($I^+$) forms
<figure><img src="image_7.png" /><figcaption>**Figure 3:** Steady-state intra-cavity pulse dynamics. Shaded areas refer to loss regions for a certain direction (striped) or both directions (checkered). (a) When the pulse is located in the middle of the cavity, gain and loss are slightly saturated and almost equal. The weak counter-propagating field experiences a significant loss. (b) At the facet of the cavity, right and left traveling intensities are equally strong. This leads to a net loss window for both directions. (c) When the pulse has left the facet, the strong gain saturation due to optical interference is still present. As the absorber is independent of this effect, a significant loss window is created.</figcaption></figure>

the pulse and grows along its path, as the right-traveling gain ($g^+$) slightly exceeds the loss everywhere. The intensity of the left-traveling wave ($I^-$) is by orders of magnitude smaller and gets further dampened as the left-traveling gain ($g^-$) is smaller than the loss. This scenario corresponds to Eq. (15) of the analytical model. When the pulse hits the facet, as in Figure 3(b), both gain components saturate to lower values than the loss, effectively creating an overall net loss window. This net loss at the facet is crucial for pulse stability and highly similar to saturable absorbing mirrors used in common passive mode-locking setups [1], [10]. At the interface, right and left traveling electric field components are present in equal strength, leading to strong optical interference. As discussed along Eq. (14), for this scenario, the saturation behavior of the gain is strongly affected by

the inversion grating term $w^+$. In contrast, the absorber in our model is only influenced by the intensity and thus saturates less. When the pulse has changed its direction, now traveling to the left, as illustrated in Figure 3(c), its tail is still prone to the overall net loss window and thus gets cut off. In this way, the pulse can sustain its shape, even in the presence of group velocity dispersion. Although this mechanism prevents the pulse from broadening, there may still be a background formation due to the finite net gain in the rest of the cavity.

This way, undesired background fluctuations may form in the middle and at the left side of the cavity, while the pulse is reshaped at the right facet. This weak field will counter-propagate toward the pulse and, at some point, reach the same location of the cavity. In this situation, the strong pulse absorbs the weak fluctuations because they experience a net loss, as described above. We emphasize that this is an intrinsic effect of the quantum system that arises when taking into account the population grating arising from optical interference in the Bloch equations. However, adding a saturable absorber with suitable properties is the crucial ingredient to generate the conditions for pulsed operation, where this effect is revealed. The described dynamics of one complete round trip are animated in Supplementary Movie 1.

Another manifestation of directional gain is symmetry-breaking in ring lasers [19], [20], [62]. There, through spontaneous emission during the formation process, one field direction randomly grows stronger, thus reducing the available gain for the other direction. With increasingly asymmetric intensities, this effect gets more pronounced, resulting in a complete decay of the counter-propagating field. Consequently, this results in unidirectional operation unless a balancing mechanism, such as backscattering, is present [20], [63].

This tendency to unidirectional operation is similar to the formation dynamics in the here-discussed case of passively mode-locked pulses, where a kind of symmetry breaking is also present. In Figure 4(a), the transient behavior of the intracavity power is depicted for the first 500 round trips over the whole cavity, starting from random noise. During the first 50 round trips, the weak initial fluctuations get slowly amplified and, due to $D_f$, filtered into a smooth but unevenly distributed intensity. Then, between ≈50 and 100 round trips, a nonlinear process takes place, rapidly amplifying randomly generated power accumulations by orders of magnitude, while hardly affecting the remaining intensity. Notably, this is only possible due to the presence of a saturable absorber, which is fast enough to follow the gain dynamics, and at the same time saturates slightly more easily. In this specific simulation, a clear pulse forms at
<figure><img src="image_8.png" /><figcaption>**Figure 4:** Visualization of (a) intracavity intensity and (b) net gain during the formation process. The zigzag trace arises from a slight intentional undersampling of the round trip frequency, such that the pulse is visualized at each point of the cavity.</figcaption></figure>

cavity positions close to 0 mm. However, a significant amount of distributed, unlocked intensity is still present in the cavity. The dynamics of background absorption and symmetry breaking are visible between round trip 100 and 200. After the background has collided several times with the pulse (twice per round trip), it has lost intensity by several orders of magnitude. In contrast to that, the pulse gains energy throughout this period of time. Therefore, these formation dynamics can be seen as a type of symmetry breaking. After the buildup process, i.e., for times >250 round trips, the pulse completely dominates the dynamics of the laser.

In Figure 4(b), the net direction-averaged gain $g_n = g_s - (a_s + a_0)$ inside the cavity is depicted. We observe similar behavior as described alongside Figure 3. In the middle of the cavity, there is a small positive net gain present, and its maximum follows the intensity peak. The loss for the counter-propagating field is not visible, as its impact on the direction-averaged gain value $g_s$ in Eq. (12) scales with its weak intensity. However, when the pulse reaches one of the facets, a significant loss window emerges, as both intensity directions significantly contribute to the saturation. As outlined above, this effect acts as a pulse shortening and stabilizing mechanism. The pulse formation dynamics sampled at each round trip are visualized in Supplementary Movie 2.

The effect of spontaneous emission noise on the laser dynamics has been tested by adding small perturbations of the electric field at every grid point, implemented similarly as in [64]. Considering the detailed active region model with

nine levels, spontaneous emission causes fluctuations of the pulse peak power with a standard deviation of ≈1 %. The observed beatnote linewidth is well below the numerical frequency resolution of ≈2 MHz. As the cavity round trip frequency is at ≈20 GHz, the timing jitter of the pulses should thus be at least 10,000 times smaller than the round trip time.

## 3.2 Investigation of operating regimes

After discussing the mode-locking dynamics in detail for a certain set of parameters, we study the robustness of the system against deviations in the absorber properties. The saturation intensities of gain and absorber ($I_{s,g}$ and $I_{s,a}$) indicate the values at which the respective variable clamps to half of its maximum. Therefore, the absorber gets saturated by a lower intracavity intensity than the gain, which is an important relation, as discussed below. Further, in the parameter sets of Tables 1 and 2, we see that the gain maximum $g_0$ is just slightly above the sum of the overall losses $a_0 + a_1$. If the difference between unsaturated gain and loss increases, the generation of passively mode-locked pulses breaks down, as background formation cannot be suppressed anymore. However, evaluating the mode-locking regime in terms of current, which is commonly used as the adjustable variable in measurements, loosens this tight condition. Assuming that in steady state, electronic transitions exactly compensate for the number of lost photons, we can calculate a photon-driven current density $J_{ph}$ in terms of the losses and intensity

$$J_{ph} = \langle (a_0 + a_s(x, t))I(x, t) \rangle \frac{e_0 L_p}{\hbar \omega_c} \quad (16)$$

Here, the overall intensity $I(x, t)$ is given as the sum of the left- and right-traveling components $I(x, t) = I^+(x, t) + I^-(x, t)$. Further, the operator $\langle \cdot \rangle$ denotes the spatial and temporal average over one round trip, and $L_p$ is the length of one QCL period. By this evaluation, we observe that a current density range of $J_{ph} \approx 15 \text{ A cm}^{-2}$ yields single pulsed states, and multipulsed states extend the narrow beatnote range by about another ≈10 A cm$^{-2}$. Therefore, even with the restriction that the unsaturated gain must not exceed the unsaturated loss too far, a significant bias-current range opens up, where the mode-locked operation may be observed experimentally. Regarding realistic QCLs, in addition to this photon-driven contribution, dark charge carrier transport holds a significant portion of the overall current. This contribution may further extend the current range where a narrow beatnote is detected, as it has only a minor influence on the field dynamics. As a comparison, the experimental device showed a threshold current density of
<figure><img src="image_9.png" /><figcaption>**Figure 5:** (a) Dependence of photon-driven current on the variation of the absorber saturation intensity. Three different regimes appear. I Single mode and weak amplitude modulation. II Single pulsed state. Decreasing the absorber intensity increases the current. Decreasing full-width-half-max (FWHM) pulse width indicates higher intensity pulses. III When the absorber bleaches too easily, unlocked multipulsing and chaotic states appear. (b) Double pulsing at the border between region II and III. (c) Chaotic emission of region III.</figcaption></figure>

≈200 A cm⁻², and a narrow beatnote regime was detected up to ≈75 A cm⁻² above this threshold. A more detailed analysis of the photon-driven current and different operational states of the simulations is given in Figure 5. Since the absorber saturation intensity *I*<sub>s,a</sub> depends on the exact overlap of the transverse optical field mode with the graphene stripes in the cavity, *I*<sub>s,a</sub> is nontrivial to determine for our one-dimensional simulation approach. Therefore, we evaluate the robustness of the system against deviations of *I*<sub>s,a</sub> while keeping the gain parameters unchanged. We obtain three different possible states that the laser may assume, and we can calculate the photon-driven current in these regimes, which is also directly proportional to the output power. In regime I, where both quantities have comparable values, the gain quickly clamps to the losses in the whole cavity, only allowing for single-mode operation or weak amplitude modulations. As this state exhibits weak power, the photon-driven current is on the order of 5 A cm⁻², and the output power some hundreds of microwatts. When the absorber bleaches ≈2.2 times more easily than the gain, the laser switches to region II, the single pulsed state. Up to a factor of 4.5 difference between *I*<sub>s,g</sub> and *I*<sub>s,a</sub> single pulse mode-locking is possible, showing that the system is tolerant

against deviations of this parameter. The associated photon-driven current covers a range between ≈10 and 25 A cm⁻². A decreasing value of the pulse's FWHM is observed when reducing *I*<sub>s,a</sub>, which goes along with higher peak intensities. Further, the pulse duration is also affected by the presence and strength of chromatic dispersion [65]–[68]. For the results shown in Figure 5(a), a value of β₂ = 5 × 10⁻²³ s² m⁻¹ was used as parameter in Eq. (1). Ideal dispersion compensation would thus lead to shorter pulses, as already mentioned in [9], while stronger dispersion broadens the pulse (compare Supplementary Information Section III.B).

Regarding the variation of saturation intensities, the case in which the absorber bleaches ≈5 times more easily than the gain results in unlocked multipulsing or chaotic behavior (region III), as depicted in Figure 5(c). However, the first data point of this region in Figure 5(a), at ≈35 A cm⁻², shows a double-pulse, see Figure 5(b), with a clear beatnote still present.

## 3.3 Transition to standard theory

We have now discussed in detail the case where the recovery times of gain and absorber are very similar, τ<sub>a</sub> ≈ τ<sub>g</sub>, as this comes closest to the experiment [9]. If the recovery time of the absorber exceeds the one of the gain (τ<sub>a</sub> ~ 10 % larger than τ<sub>g</sub>), the rapid nonlinear amplification, visible in Figure 4(a) between ≈50 and 100 round trips, cannot take place and, therefore, pulsed operation is not supported. We note that the ratio of recovery times may change due to different possible realizations of gain media (THz, mid-IR QCLs, IC, or QD lasers) and possible improvements in the recovery speed of saturable absorbers. Furthermore, different values of gain recovery times in THz QCLs have been reported in literature, ranging from 2 ps to several tens of ps [15], [69]–[71]. In the Supplementary Information, we depict the maximum ratio of gain and loss over the ratio of gain and absorber recovery times that still yields stable pulses. In order to extend the coverage of our model for the cases with longer gain recovery time, we show results of pulsed operation for τ<sub>g</sub> = 10.5 ps and τ<sub>a</sub> = 3 ps in Figure 6(a). The remaining used parameters, as well as a color map of the spatiotemporal intensity, can be found in the Supplementary Information Section III.C. Notably, in this situation, there is a large loss window behind the pulse, allowing for very strong background absorption. Thus, a region with significant gain in front of the pulse can be balanced without yielding chaotic behavior. In fact, this specific simulation was carried out using a₀ = 11.0 cm⁻¹, such that the gain maximum exceeds the unsaturated loss by g₀ − (a₀ + a₁) = 1 cm⁻¹. Therefore, assuming this rather long gain recovery time, the model reacts less
<figure><img src="image_10.png" /><figcaption>**Figure 6:** Converged pulse when the absorber recovers significantly faster than the gain. (a) The gain recovers within 10.5 ps, the absorber within 3 ps. (b) The gain recovers within 2.6 ps, the absorber within 1 ps. In both cases, a net loss window is always present behind the pulse.</figcaption></figure>

sensitively to the gain/loss ratio than for similar recovery times, as described above. This observation goes along with the discussion at the very beginning of Section 2, as this simulation rather approaches the scenario of intermediate gain recovery times. Hence, for $\tau_g = 10.5$ ps, the strong demands on the gain and absorber parameters may be significantly relaxed, while pulses can still be obtained.

Additionally, there is the possibility left where the absorber is significantly faster than the gain, $\tau_a < \tau_g$, while the gain recovery time remains at $\tau_g = 2.6$ ps. Therefore, we set $\tau_a = 1$ ps and leave the remaining parameters as in Tables 1 and 2. In Figure 6(b), the resulting steady-state distributions of intensity and gain are depicted. A single pulse emerges that is significantly shorter, with higher peak intensity, as compared to the initial simulation (Figure 3). In Figure 6(b), a net loss window emerges behind the pulse. Therefore, the tail of the pulse gets constantly cut off, while the leading edge sees a large amount of gain. Thus, photons are created at the front of the pulse and get re-absorbed at its end. As a result of this asymmetry, the pulse velocity is higher than the intrinsic group velocity of the medium, which determines the speed of background fluctuations. The dominating dynamics of this case are highly similar to the one of laser self-pulsations, described in [72]. The superluminal motion of the pulse additionally enables the loss window to suppress the amplification of background

fluctuations throughout the whole cavity and, thus, to sustain pulsed operation.

## 4 Summary and conclusion

To summarize, we have presented a model of single pulse passive mode-locking in lasers with fast recovery time, i.e., faster than the cavity round trip time, such as THz QCLs. The presence of an incoherent distributed saturable absorber, realized by graphene stripes embedded along the top contact, can lead to pulse formation if it saturates with a suitable intensity and recovers at least equally fast as the gain. The incoherent optical properties of graphene and its high charge carrier mobility, suppressing SHB, make it an ideal choice for the absorber and enable the formation of a net loss window at the cavity facets. An optical interference term appearing in the gain description, associated with SHB, enables the suppression of a continuous wave background. The model shows a certain robustness against parameter variations and allows mode-locking in a significant current density bias range. Furthermore, the case when the absorber is faster than the gain is investigated. A shorter pulse emerges, which is followed by a pronounced loss window. Thus, the stabilizing effect is not restricted to the cavity interfaces, and interference effects become less important. Our work increases the understanding of experimentally realized passive pulse formation in lasers with fast gain dynamics, such as THz QCLs, and provides a guide on how shorter pulses can be achieved. The described model may also be relevant for the design of other mode-locked lasers with quantum-engineered active regions, such as IC and QD lasers.

**Research funding:** The authors acknowledge financial support by the European Union's QuantERA II [G.A. n. 101017733] – QATACOMB Project “Quantum correlations in terahertz QCL combs” (Funding organizations: DFG – Germany [Project n. 491801597]), European Union's Horizon 2020 research and innovation program under grant agreement No. 964735 (FET-OPEN EXTREME-IR).

**Author contributions:** All authors have accepted responsibility for the entire content of this manuscript and approved its submission.

**Conflict of interest:** Authors state no conflicts of interest.

**Data availability:** The datasets generated and/or analyzed during the current study are available from the corresponding author upon reasonable request.
References

[1] H. Haus, "Mode-locking of lasers," *IEEE J. Sel. Top. Quantum Electron.*, vol. 6, no. 6, pp. 1173-1185, 2000.

[2] T. Udem, R. Holzwarth, and T. W. Hänsch, "Optical frequency metrology," *Nature*, vol. 416, no. 6877, pp. 233–237, 2002.

[3] T. W. Hänsch, "Nobel lecture: passion for precision," *Rev. Mod. Phys.*, vol. 78, no. 4, pp. 1297–1309, 2006.

[4] A. H. Zewail, "Femtochemistry: atomic-scale dynamics of the chemical bond," *J. Phys. Chem. A*, vol. 104, no. 24, pp. 5660–5694, 2000.

[5] S. A. Diddams, K. Vahala, and T. Udem, "Optical frequency combs: coherently uniting the electromagnetic spectrum," *Science*, vol. 369, no. 6501, p. eay3676, 2020.

[6] B. Schwarz, et al., "Monolithically integrated mid-infrared lab-on-a-chip using plasmonics and quantum cascade structures," *Nat. Commun.*, vol. 5, no. 1, p. 4085, 2014.

[7] S. Bartalini, et al., "Frequency-comb-assisted terahertz quantum cascade laser spectroscopy," *Phys. Rev. X*, vol. 4, no. 2, p. 021006, 2014.

[8] T. Fortier and E. Baumann, "20 years of developments in optical frequency comb technology and applications," *Commun. Phys.*, vol. 2, no. 1, p. 153, 2019.

[9] E. Riccardi, et al., "Short pulse generation from a graphene-coupled passively mode-locked terahertz laser," *Nat. Photon*, vol. 17, no. 7, pp. 607–614, 2023.

[10] U. Keller, et al., "Semiconductor saturable absorber mirrors (SESAM's) for femtosecond to nanosecond pulse generation in solid-state lasers," *IEEE J. Sel. Top. Quantum Electron.*, vol. 2, no. 3, pp. 435–453, 1996.

[11] H. W. Mocker and R. J. Collins, "Mode competition and self-locking effects in a Q-switched ruby laser," *Appl. Phys. Lett.*, vol. 7, no. 10, pp. 270–273, 1965.

[12] A. Komarov, H. Leblond, and F. Sanchez, "Passive harmonic mode-locking in a fiber laser with nonlinear polarization rotation," *Opt. Commun.*, vol. 267, no. 1, pp. 162–169, 2006.

[13] H. A. Haus, J. G. Fujimoto, and E. P. Ippen, "Structures for additive pulse mode locking," *JOSA B*, vol. 8, no. 10, pp. 2068–2076, 1991.

[14] J. Faist, F. Capasso, D. L. Sivco, C. Sirtori, A. L. Hutchinson, and A. Y. Cho, "Quantum cascade laser," *Science*, vol. 264, no. 5158, pp. 553–556, 1994.

[15] J. Riepl, et al., "Field-resolved high-order sub-cycle nonlinearities in a terahertz semiconductor laser," *Light Sci. Appl.*, vol. 10, no. 1, p. 246, 2021.

[16] A. Hugi, G. Villares, S. Blaser, H. C. Liu, and J. Faist, "Mid-infrared frequency comb based on a quantum cascade laser," *Nature*, vol. 492, no. 7428, pp. 229–233, 2012.

[17] M. Piccardo, et al., "Frequency combs induced by phase turbulence," *Nature*, vol. 582, no. 7812, pp. 360–364, 2020.

[18] M. Piccardo and F. Capasso, "Laser frequency combs with fast gain recovery: physics and applications," *Laser Photonics Rev*, vol. 16, no. 2, p. 2100403, 2022.

[19] B. Meng, M. Singleton, J. Hillbrand, M. Franckié, M. Beck, and J. Faist, "Dissipative Kerr solitons in semiconductor ring lasers," *Nat. Photon*, vol. 16, no. 2, pp. 142–147, 2022.

[20] L. Seitner, et al., "Backscattering-induced dissipative solitons in ring quantum cascade lasers," *Phys. Rev. Lett.*, vol. 132, no. 4, p. 043805, 2024.

[21] Q. Lu, F. Wang, D. Wu, S. Slivken, and M. Razeghi, "Room temperature terahertz semiconductor frequency comb," *Nat. Commun.*, vol. 10, no. 1, p. 2403, 2019.

[22] L. Consolino, et al., "Direct observation of terahertz frequency comb generation in difference-frequency quantum cascade lasers," *Appl. Sci.*, vol. 11, no. 4, p. 1416, 2021.

[23] J. Popp, et al., "Self-consistent simulations of intracavity terahertz comb difference frequency generation by mid-infrared quantum cascade lasers," *J. Appl. Phys.*, vol. 133, no. 23, 2023, https://doi.org/10.1063/5.0151036.

[24] D. Burghoff, et al., "Terahertz laser frequency combs," *Nat. Photon*, vol. 8, no. 6, pp. 462–467, 2014.

[25] C. R. Menyuk and M. A. Talukder, "Self-induced transparency modelocking of quantum cascade lasers," *Phys. Rev. Lett.*, vol. 102, no. 2, p. 023903, 2009.

[26] P. Tzenov, et al., "Passive and hybrid mode locking in multi-section terahertz quantum cascade lasers," *New J. Phys.*, vol. 20, no. 5, p. 053055, 2018.

[27] M. A. Talukder and C. R. Menyuk, "Quantum coherent saturable absorption for mid-infrared ultra-short pulses," *Opt. Express*, vol. 22, no. 13, pp. 15608–15617, 2014.

[28] A. M. Perego, B. Garbin, F. Gustave, S. Barland, F. Prati, and G. J. De Valcárcel, "Coherent master equation for laser modelocking," *Nat. Commun.*, vol. 11, no. 1, p. 311, 2020.

[29] J. Hausen, K. Lüdge, S. V. Gurevich, and J. Javaloyes, "How carrier memory enters the Haus master equation of mode-locking," *Opt. Lett.*, vol. 45, no. 22, pp. 6210–6213, 2020.

[30] L. Columbo, et al., "Unifying frequency combs in active and passive cavities: temporal solitons in externally driven ring lasers," *Phys. Rev. Lett.*, vol. 126, no. 17, p. 173903, 2021.

[31] F. Bianco, et al., "THz saturable absorption in turbostratic multilayer graphene on silicon carbide," *Opt. Express*, vol. 23, no. 9, pp. 11632–11640, 2015.

[32] V. Bianchi, et al., "Terahertz saturable absorbers from liquid phase exfoliation of graphite," *Nat. Commun.*, vol. 8, no. 1, p. 15763, 2017.

[33] C. Jirauschek, M. Riesch, and P. Tzenov, "Optoelectronic device simulations based on macroscopic Maxwell-Bloch equations," *Adv. Theory Simul.*, vol. 2, no. 8, p. 1900018, 2019.

[34] H. A. Haus, "Theory of mode locking with a fast saturable absorber," *J. Appl. Phys.*, vol. 46, no. 7, pp. 3049–3058, 1975.

[35] F. X. Kärtner, J. Aus Der Au, and U. Keller, "Mode-locking with slow and fast saturable absorbers — what's the difference?" *IEEE J. Sel. Top. Quantum Electron.*, vol. 4, no. 2, pp. 159–168, 1998.

[36] N. Opačak and B. Schwarz, "Theory of frequency-modulated combs in lasers with spatial hole burning, dispersion, and Kerr nonlinearity," *Phys. Rev. Lett.*, vol. 123, no. 24, p. 243902, 2019.

[37] N. Opačak, S. D. Cin, J. Hillbrand, and B. Schwarz, "Frequency comb generation by Bloch gain induced giant Kerr nonlinearity," *Phys. Rev. Lett.*, vol. 127, no. 9, p. 093902, 2021.

[38] D. Burghoff, "Unraveling the origin of frequency modulated combs using active cavity mean-field theory," *Optica*, vol. 7, no. 12, pp. 1781–1787, 2020.

[39] C. Y. Wang, et al., "Coherent instabilities in a semiconductor laser with fast gain recovery," *Phys. Rev. A*, vol. 75, no. 3, p. 031802(R), 2007.

[40] K. I. Bolotin, et al., "Ultrahigh electron mobility in suspended graphene," *Solid State Commun.*, vol. 146, nos. 9–10, pp. 351–355, 2008.
[41] L. Banszerus, et al., "Ballistic transport exceeding 28 µm in CVD grown graphene," Nano Lett., vol. 16, no. 2, pp. 1387-1391, 2016.

[42] D. R. Cooper, et al., "Experimental review of graphene," *ISRN Condens. Matter Phys.*, vol. 2012, pp. 1–56, 2012.

[43] F. Bonaccorso, Z. Sun, T. Hasan, and A. Ferrari, "Graphene photonics and optoelectronics," *Nat. Photon.*, vol. 4, no. 9, pp. 611–622, 2010.

[44] H. A. Hafez, S. Kovalev, K.-J. Tielrooij, M. Bonn, M. Gensch, and D. Turchinovich, "Terahertz nonlinear optics of graphene: from saturable absorption to high-harmonics generation," *Adv. Optical Mater.*, vol. 8, no. 3, p. 1900771, 2020.

[45] Q. Bao, et al., "Atomic-layer graphene as a saturable absorber for ultrafast pulsed lasers," *Adv. Funct. Mater.*, vol. 19, no. 19, pp. 3077-3083, 2009.

[46] G. Xing, H. Guo, X. Zhang, T. C. Sum, and C. H. A. Huan, "The physics of ultrafast saturable absorption in graphene," *Opt. Express*, vol. 18, no. 5, pp. 4564-4573, 2010.

[47] C. Jirauschek, "Theory of hybrid microwave-photonic quantum devices," *Laser Photonics Rev.*, vol. 17, no. 12, p. 2300461, 2023.

[48] L. Allen and J. H. Eberly, *Optical Resonance and Two-Level Atoms*, New York, Dover Publications, 1987.

[49] C. Jirauschek and T. Kubis, "Modeling techniques for quantum cascade lasers," *Appl. Phys. Rev.*, vol. 1, no. 1, p. 011307, 2014.

[50] L. A. Sterczewski, et al., "Interband cascade laser frequency combs," *J. Phys. Photonics*, vol. 3, no. 4, p. 042003, 2021.

[51] P. Bardella, L. L. Columbo, and M. Gioannini, "Self-generation of optical frequency comb in single section quantum dot Fabry-Perot lasers: a theoretical study," *Opt. Express*, vol. 25, no. 21, pp. 26234-26252, 2017.

[52] L. Lugiato, F. Prati, and M. Brambilla, *Nonlinear Optical Systems*, Cambridge, Cambridge University Press, 2015.

[53] R. W. Boyd, *Nonlinear Optics*, 4th ed. London, Academic Press, 2020.

[54] T. S. Mansuripur, et al., "Single-mode instability in standing-wave lasers: the quantum cascade laser as a self-pumped parametric oscillator," *Phys. Rev. A*, vol. 94, no. 6, p. 063807, 2016.

[55] M. Piccardo, et al., "Time-dependent population inversion gratings in laser frequency combs," *Optica*, vol. 5, no. 4, pp. 475-478, 2018.

[56] S. Meinecke, L. Drzewietzki, C. Weber, B. Lingnau, S. Breuer, and K. Lüdge, "Ultra-short pulse generation in a three section tapered passively mode-locked quantum-dot semiconductor laser," *Sci. Rep.*, vol. 9, no. 1, p. 1783, 2019.

[57] D. Brida, et al., "Ultrafast collinear scattering and carrier multiplication in graphene," *Nat. Commun.*, vol. 4, no. 1, p. 1987, 2013.

[58] C. Jirauschek, "Density matrix Monte Carlo modeling of quantum cascade lasers," *J. Appl. Phys.*, vol. 122, no. 13, p. 133105, 2017.

[59] H. Li, et al., "Graphene-coupled terahertz semiconductor lasers for enhanced passive frequency comb operation," *Adv. Sci.*, vol. 6, no. 20, p. 1900460, 2019.

[60] F. P. Mezzapesa, et al., "Terahertz frequency combs exploiting an on-chip, solution-processed, graphene-quantum cascade laser coupled-cavity," *ACS Photonics*, vol. 7, no. 12, pp. 3489-3498, 2020.

[61] B. S. Williams, "Terahertz quantum-cascade lasers," *Nat. Photon.*, vol. 1, no. 9, pp. 517-525, 2007.

[62] E. J. D'Angelo, E. Izaguirre, G. B. Mindlin, G. Huyet, L. Gil, and J. R. Tredicic, "Spatiotemporal dynamics of lasers in the presence of an imperfect O(2) symmetry," *Phys. Rev. Lett.*, vol. 68, no. 68, pp. 3702-3705, 1992.

[63] N. Opačak, et al., "Nozaki-Bekki solitons in semiconductor lasers," *Nature*, vol. 625, no. 7996, pp. 685-690, 2024.

[64] L. L. Columbo, P. Bardella, and M. Gioannini, "Self-pulsing in single section ring lasers based on quantum dot materials: theory and simulations," *Opt. Express*, vol. 26, no. 15, pp. 19044-19058, 2018.

[65] F. Wang, et al., "Short terahertz pulse generation from a dispersion compensated modelocked semiconductor laser," *Laser Photonics Rev.*, vol. 11, no. 4, p. 1700013, 2017.

[66] M. A. Talukder and C. R. Menyuk, "Self-induced transparency modelocking of quantum cascade lasers in the presence of saturable nonlinearity and group velocity dispersion," *Opt. Express*, vol. 18, no. 6, pp. 5639-5653, 2010.

[67] P. Tzenov, D. Burghoff, Q. Hu, and C. Jirauschek, "Time domain modeling of terahertz quantum cascade lasers for frequency comb generation," *Opt. Express*, vol. 24, no. 20, pp. 23232-23247, 2016.

[68] L. Seitner, J. Popp, M. Riesch, M. Haider, and C. Jirauschek, "Group velocity dispersion in terahertz frequency combs within a generalized Maxwell-Bloch framework," *J. Phys.: Conf. Ser.*, vol. 2090, no. 1, p. 012082, 2021.

[69] R. P. Green, et al., "Gain recovery dynamics of a terahertz quantum cascade laser," *Phys. Rev. B*, vol. 80, no. 7, p. 075303, 2009.

[70] D. R. Bacon, et al., "Gain recovery time in a terahertz quantum cascade laser," *Appl. Phys. Lett.*, vol. 108, no. 8, p. 081104, 2016.

[71] C. G. Derntl, et al., "Gain dynamics in a heterogeneous terahertz quantum cascade laser," *Appl. Phys. Lett.*, vol. 113, no. 18, p. 181102, 2018.

[72] H. Risken and K. Nummedal, "Self-pulsing in lasers," *J. Appl. Phys.*, vol. 39, no. 10, pp. 4662-4672, 1968.

**Supplementary Material:** This article contains supplementary material (https://doi.org/10.1515/nanoph-2023-0657).