1. Introduction

Electrons provide the fundamental first step in response of matter to light. The feasibility of shaping light pulses at the scale of individual oscillations, from mid-IR to UV [1,2], offers rich opportunities for controlling electronic response to light on sub-cycle timescale (e.g. [3–13]), leading to a variety of fascinating phenomena such as optically induced anomalous Hall effect [14–16], topological phase transitions with polarization-tailored light [13], or the topological resonance [7]. Over multiple laser cycles, control of electron dynamics with light also enables the so-called Floquet engineering – the tailored modification of the cycle-average properties of a light-dressed system, see e.g. [17] for a recent review.

In this context, starting with the pioneering work [18], high harmonic spectroscopy has developed into a powerful tool for exploring laser-driven electron dynamics in solids, see e.g. recent reviews [19–21]. Examples include identification of the common physical mechanisms underlying high harmonic generation in atoms, molecules and solids (e.g. [10,22]), observation of Bloch oscillations [23], resolving interfering pathways of electrons in crystals with about 1-fsec precision [24], inducing [13] and monitoring topological [25–27] and Mott insulator-to-metal [28] phase transitions, resolving coherent oscillation of electronic charge density order [29], identifying the van Hove singularities in the conduction bands [30], and reconstructing effective potentials seen by the light-driven electrons with picometer accuracy [9].

Here we apply machine learning to the analysis of high harmonic generation from a crystal, which allows us to reconstruct the band structure of the crystal and characterize incident few-cycle laser pulses, including both their nonlinear chirp and the phase of the carrier oscillations under the envelope (CEP).

The fundamental role of the CEP in nonlinear light-matter interaction has been understood theoretically in [31–35], stimulating first experiments in the gas phase [36]. Gas-phase methods for characterizing few-cycle pulses and their CEP are generally based on photo-electron spectroscopy and include stereo-above-threshold ionization (stereo-ATI) [37–39] and attosecond streak camera with its modifications [40–44]. However, these methods are complex and require sophisticated vacuum-based measurement techniques, which have to be further added to the experimental setup. As a result, such measurements are currently only available in a very few laboratories world-wide. Using nonlinear response of solids for characterizing the CEP has also been pursued [45–47]. The pioneering experiment [48] that has been just reported shows the importance of developing all-solid state recognition systems which characterize the pulse and also the material in one go. This work shows that the waveform can be extracted using nonlinear response, subject to well-known and well-characterized nonlinearity, for relatively weak pulses, and represents a breakthrough in the development of all solid-state pulse characterization methods.

The situation becomes more complex in the regime of intense light-matter interaction, where the nonlinearity is often unknown. We extend the analysis to this high-intensity regime, looking at the pulse directly in the interaction region, and assuming no detailed apriori knowledge of the material. We therefore suggest a route to developing a solid-state based method which enables complete pulse characterization and simultaneously gives access to the effective Hamiltonian of the crystal. By implication, the also gives access to the dynamics underlying the nonlinear optical response. Particularly relevant to our work is the use of two-color high harmonic spectroscopy for all-optical reconstruction of the band structure [49] from the two-dimensional high harmonic spectra, recorded as a function of the harmonic frequency and the two-color delay.

Two-dimensional spectra of the highly nonlinear-optical response often provide sufficient information to recover the pulse. If such recovery is placed on a solid footing, it would allow one to characterize the laser pulse directly in the interaction region – a sought-after goal. In special cases, such as the case of the two-color high harmonic spectroscopy of attosecond pulses using fundamental and the second harmonic (e.g. [50]), the 2D harmonic spectra recorded as a function of the two-color phase and the harmonic frequency does carry sufficient information for reconstructing attosecond pulses. Such reconstruction does, however, require detailed understanding of the physics of the microscopic quantum response. The possibility of solving a full reconstruction problem in a general case, recovering both the pulse and the quantum system, remains unexplored.

Application of neural networks is a natural way to address this problem. Indeed, when there is no physically transparent, simple, or well known functional dependence between the response data and the parameters one wishes to reconstruct, the problem is well suited for neural networks. Such networks aim to find a smooth analytical function $f_\theta (x)$ which connects the input $x_i$ and the desired output $y_i$. If it is successful, one can conclude that the data $x$ indeed does contain the information $y$. Pertinent examples include applications to solving the Schrödinger equation, where neural networks can output highly accurate results [51,52]. However, the previous approaches do not aim to reconstruct the dynamical properties of quantum systems from macroscopic measurements. Our results show that, given sufficient training set, the 2D spectra of the high harmonic response as a function of the nonlinear response frequency and the CEP of an unknown driving pulse allow for simultaneous reconstruction of both the pulse and the unknown crystal band structure, defined by a standard set of parameters such as on-site energies and hoppings, even in the presence of substantial noise including CEP jitter.

We stress that our goal here is not to develop new machine learning methods, but rather to apply the already established methods to challenging problems in attosecond physics of solids and extreme nonlinear optics. Thus, we do not focus on the detailed analysis of pro-s and con-s of particular neural networks in favour of standard tools.

2. Problem description

To demonstrate the method, we assume that the nonlinear medium can be described by the tight-binding approximation (TBA) up to the 6th hopping order, and the polynomial expansion of pulse’s spectral phase around $\omega _0$ includes terms up to the 3rd (i.e. the CEP, the chirp, and the cubic phase; the linear term which amounts to a simple pulse delay is dropped). Neural networks are trained to recover these parameters. Recovery of coefficients for the TBA is chosen instead of direct recovery of the band shape due to TBA’s robustness and flexibility, used in e.g. the Wannierization procedure implemented in the DFT code "Quantum Espresso" [53]. This makes our approach general, flexible, and transferable to many systems. The same philosophy is applied to recovering the parameters of the driving laser pulse.

For the quantum system, we begin with a modified Rice-Mele model [54] with nearest neighbor, next nearest neighbor, etc. hoppings, see Fig. 1(a) (and Supplement 1 information for further details). Both the on-site energies and the couplings are assumed to be unknown. The reconstruction procedure is expected to output both the parameters of the pulse and of the lattice. We are thus faced with a nonlinear optimization problem with a very large input dimension, which requires finding optimal interpolation between the existing trial samples.

 

Fig. 1. Model system used for pulse reconstruction and band structure spectroscopy via nonlinear-optical response. (a) The two types of sites present, A and B, are connected by hopping constants $t_j$ between similar sites and $h_1$ between the next-neighbor sites of the different type. (b-d) Examples of band structures that can be generated with this model system.

Download Full Size | PDF

Such a regression tool is provided by deep neural networks (DNNs) [55], already used for such diverse applications as boosting the signal-to-noise ratio in LHC collision data [56], establishing a fast mapping between galaxy and dark matter distribution [57], and constructing efficient representations of many-body quantum states [58]. In a recent work [59], deep neural networks were used to classify molecules and recover their parameters using ultrafast dipole responses in the time domain. This work suggests that neural networks can also be trained to classify molecules and reconstruct their parameters using HHG intensity spectra, which is more compatible with experiment.

The inherent resilience of neural networks to noise is an important asset for pulse shape characterization. The emergence of photonic implementations of feed-forward neural networks [60] outlines a perspective of implementing this regression scheme in an all-optical way.

In our simulations, the vector potential of the incident laser field, $A(t)$, is generated in the frequency domain with the unknown to the neural network quadratic and the cubic phases for the carrier frequency $\omega _0$ and the pulse width $\sigma =2\pi /\omega _0$:

Three different datasets are used for our reconstruction procedure. The first one contains pulses with a flat spectral phase ($\beta =\epsilon =0$) and varying CEP, the second with varying chirp and CEP ($\epsilon =0$), while for the third dataset we vary all parameters of the pulse as outlined above. The typical pulses present in the respective datasets are shown in Fig. 2, both in time and frequency domain. The latter shows the spectral phases (red curves) alongside the spectral amplitudes (blue curves), as a function of $\omega /\omega _0$, where $\omega _0$ is the carrier. Our typical simulated "measurement" assumes that one can systematically vary the (unknown) CEP. Thus, we perform calculations by varying $\varphi + \Delta \varphi$, with $\Delta \varphi$ spanning the full range $\Delta \varphi \in \left [0,2\pi \right )$ in $N_\varphi =32$ steps. For each $\Delta \varphi$, we measure the absolute value of the spectral amplitude of the laser-induced current $|j(\omega )|$, which is given by the Fourier transform of the calculated current $j(t)$.

 

Fig. 2. Sample pulse shapes used in the problem. (a-c) show the time-domain. (d-f) are the frequency-domain representations of their respective pulses, where the spectral phase (red line) and spectral amplitude normalized to the central frequency amplitude (blue area), are plotted with respect to frequency. The pulse parameters are $(\varphi, \beta, \epsilon )$, left to right: $(\pi /2, 2.0, 0.0), (0.0, 2.0, -0.5), (0.0, 1.0, 1.0)$.

Download Full Size | PDF

The input data is composed of the absolute values of the integer harmonic amplitudes $|j(N\omega, \varphi +\Delta \varphi )|$. The resulting 2D map as a function of $\Delta \varphi$ and $\omega$ is used as the input into the neural network. The network must then infer the (randomly chosen in each trial) intraband hoppings $\{t_j\}$, the unknown initial CEP $\varphi$, and the pulse parameters $\alpha$ and $\lambda$. The values of $\sigma$, frequency $\omega$, and $h_1$ are kept constant throughout. For more information about the inputs used, see Supplement 1 information.

3. Methods

In the reconstruction procedure, we have used two separate neural networks. The first was used for recovering those parameters that should remain invariant as the data is changed along the $\varphi$ axis: these are the band and spectral phase expansion coefficients. The second network was used for recovering non-invariant parameters, i.e. the CEP. Both are constructed using the same architecture.

The building blocks of this architecture consist of two alternating transformations performed upon the 4-index tensors $x_{klmn}$. Here the index $k$ labels the frequency, $l$ labels the CEP angle, $m$ is the channel number (there’s only 1 channel in the input data), and $n$ is the sample number.

The first kind of transformation, which is depicted as vertical purple stripes in Fig. 3, is a vertical 1-dimensional convolution along $\omega$. The reduction of resolution along the $\omega$ axis with increasing layer depth (downsampling, aka pooling; see Fig. 3) is effected by strided convolutions (see Supplement 1for further details).

 

Fig. 3. The network architecture used in the work. For the CEP resolution problem, the symmetry layer has complex outputs, and the final FC layer is complex-valued. The numbers under the blocks denote the number of channels in the respective layers; the tilted numbers show the size of the 1D convolutional filters.

Download Full Size | PDF

The second kind of transformation, called the circulant layer, is depicted as horizontal brown stripes in Fig. 3. It multiplies the data $\mathbf {x}$ by a 4-index tensor $\mathbf {C}$. Each slice of the C-tensor $\mathbf {C}_{\cdot \cdot k l}$ is a circulant matrix which commutes with any circular shifts along $\varphi$, thereby preserving the invariance of the data to such shifts. The transformation has the following form:

Each of these linear transformations is preceded by a batch normalization [61] layer, and the Mish activation [62] is applied after the circulant layers, see Fig. 3.

For the first neural network, which should recover the parameters that are not sensitive to the CEP of the laser pulse, we force the required invariance of the output by summing the outputs of the last circulant layer along the $\varphi$ dimension (2nd index) before applying the dense layers, therefore rendering the exact positions of the values unavailable to the NN:

The output data $y_{km;n}$ is then acted upon by a linear transform defined by a real-valued matrix $W_{j;km}$ before being compared to the correct shift-invariant parameters $p_{jn}$.

For the second neural network, we apply a modified symmetry layer to encode the position of the abstract features defined by the initial layers of the NN in the phase of the first harmonic of FFT:

Using circulant matrices and symmetry layers therefore allows one to introduce circular shift invariance into the problem explicitly, rather than rely on the neural network learning it by itself. By itself, this could be achieved by using convolutional layers and circular padding. However, the procedure we introduce allows one to construct mappings with an infinite "field of vision" where the output of a particular pixel $x_{klmn}$ depends on the inputs at each phase angle $\varphi$, as opposed to some finite window $[\varphi _l-2\pi N_l/N_\varphi ; \varphi _l+2\pi N_l/N_\varphi ]$ as is the case for convolutional networks, where $N_l$ is the convolutional filter size along $l$ (see Supplement 1). Indeed, one could expect, e.g. from semiclassical considerations, that the closer two frequencies are, the more similar are the electron dynamics behind their emission, and therefore the correlations found between signals at neighboring frequencies are more likely to be meaningful than those between distant ones. However, the same cannot be said about the CEP angle.

For either problem, this network is trained for 200 epochs with the AdaBelief [63] optimizer with batch size 256; the learning rate is initially set to $10^{-3}$ and discounted at the 190th epoch by a factor of 10.

For the Rice-Mele model, we generate three distinct datasets to evaluate the performance of the neural network in reconstructing pulses with varying complexity and varying influence of the CEP on the nonlinear response. Each dataset includes 262144 samples, each sample consisting of a 41x32 HHG spectrum, and the correct pulse and band parameters. The first dataset contains pulses with an unknown CEP and a flat spectral phase and crystals with 6 unknown band parameters. The second one contains pulses with a random chirp and 5 band parameters, the third contains pulses with a random chirp and a cubic phase, and 4 band parameters. The number of band parameters is varied to keep the dataset density the same for a more informative comparison, see Supplement 1, Section 7, for more information on the datasets used. The reconstruction performance is shown in tables S1-S3, respectively.

To speed up training, we train networks to work on more complex datasets (e.g. with chirped pulses) by using a network pre-trained on a simpler dataset (e.g. with no chirp).

To simulate experimental noise, we have introduced a random CEP shift to each of the pulses in each sample, sampled from a uniform distribution within $-50\div 50\text { mrad}$. For each sample, the output amplitudes were also affected by a uniformly-distributed multiplicative noise with an amplitude of $1\%$ times the maximal amplitude of the 3rd harmonic registered for any CEP for the given sample; this allowed to simulate a spectrometer with a dynamic range of 2 orders of magnitude.

4. Results

The recovery results for the first dataset are demonstrated in Fig. 4. The agreement between the actual and the reconstructed data, both for the system and for the pulse, is excellent. We observe, however, that it is not uniform across the entire range of parameters. The ability of the neural network to reconstruct pulse and material parameters benefits from the dynamic range and the spectral width of the recorded signal. Conversely, as discussed in Supplement 1, the effective shrinking of the HHG spectrum caused by background noise and large spectral phase degrades its accuracy.

 

Fig. 4. Band parameters (a-f) and CEP (g, h) recovered by the neural network for the source problem (Rice-Mele model, 1st dataset, no chirp or cubic phase). The solid orange line $y=x$ serves as a reference. All points used belong to the test set.

Download Full Size | PDF

After demonstrating that our neural network recovers the band structure and pulse parameters of the simple source (Rice-Mele) model with high accuracy, we applied our approach to a more complex system. As an example, we used the 2-dimensional gapped graphene system, and generated another smaller dataset organized in the same way as before. Its parameters (on-site energy and first-neighbor hopping) were generated in the vicinity of the experimental values for hBN, namely within $4.0\div 7.27$ eV and $0.08\div 0.16$ a.u., respectively. A more detailed description of the computational methods used for this model can be found in Supplement 1. We simulated its response by integrating the semiconductor Bloch equations. To simulate experimental conditions, we applied noise using the same procedure as described above.

We applied the transfer learning approach by using the networks pre-trained on the dataset with no quadratic or cubic phase, and partially retraining them (see Supplement 1). We discovered that, in spite of the greater complexity of the underlying physical system, such a setup recovers the CEP with a precision comparable to the original problem (see Figs. 5), and achieves good relative accuracy on the band parameters.

 

Fig. 5. Band structure (a), pulse waveform (b), and CEP (c, d) recovered by the neural network for the target problem (2-dimensional gapped graphene with phase and multiplicative noise).

Download Full Size | PDF

This adds new significance to the achieved results. Indeed, we have now demonstrated that the spectra from the initial model, while only resembling a real-world setup in a qualitative sense, provide a useful training ground for the neural network, allowing it to acquire useful abstract concepts (parametrized by the deeper layers) which it can later apply to more practical problems, for which it was also harder to generate as many training samples.

5. Discussion

It is hardly surprising that solid-state HHG spectra contain information on both the pulse parameters and the parameters of the quantum system. Our results show, however, that this information is sufficiently rich to enable reconstruction of both the laser pulse parameters, such as spectral phase, and the crystal parameters in terms of its parametrized band structure. We found that parametrization of both the laser pulse and the crystal structure were crucial for robust reconstruction, reducing the number of unknown parameters.

Our approach replicates the advantages of gas-phase HHG spectroscopy, namely, the ability to resolve the CEP of laser pulses and the complex pulse shapes with polynomial spectral phase. At the same time, it requires neither the XUV pulses nor the photoelectron spectroscopy, such as the stereo-ATI. Moreover, it allows for all-optical solid-state implementation. In terms of required observables, it is closest to the method based on measuring the half-cycle cutoffs in gas-phase high harmonic generation [64]. However, it also allows one to deal with very strong chirps. Applying neural networks to the analysis of half-cycle cutoffs and high harmonic generation spectra in the gas phase could be very interesting, especially in molecules, where multiple coupled harmonic generation channels present challenges for unravelling the underlying laser-driven multi-electron dynamics [65].

One could apply the developed neural network design to standard gas-phase experiments, to analyse the possibility of resolving the spectral phase and the CEP for short pulses by processing HHG spectra generated by known inert gases (Ar, Ne, etc.). In this case the neural network can be trained using TDSE simulations of the necessary responses before being applied to real experimental data.

The key difficulty of using high harmonic spectroscopy in solids is that, without apriori knowledge of the band structure, one lacks closed-form solutions for electron dynamics, similar to those available in the gas-phase. Our method circumvents this difficulty. Pulse characterization device implementing our principles could be tabletop, all solid-state, and capable of operating at ambient conditions.

Another interesting direction to pursue would be to apply novel physics-informed neural network architectures [66,67] to resolve Hamiltonians of systems with many degrees of freedom (such as molecules) using time-resolved HHG spectra, such as those obtained from solids driven by mid-IR fields [24]. Neural networks can also be used for processing the sets of harmonic spectra connected by other relations, such as being measured for different angles between the crystal axes and the driving field, to uncover effective laser-modified potentials for the charge motion, extending the pioneering work in Ref. [9] to recover effective potentials of active band electrons.

The obvious next step of this work is to validate our proposal experimentally. It is important to outline what such experimental validation needs. Any proposed new method must be calibrated against well established ones. In our case, this means performing dedicated experiments on solid-state HHG with the simultaneous measurement of the driving pulses and the emitted XUV light with an attosecond streak-camera. The experiment must use CEP-controlled incident few-cycle pulse and perform systematic and extensive measurements of the HHG response vs. CEP while varying quadratic and cubic chirps in controlled manner, for a range of incident pulse energies. While there are currently no such experimental data available, the required setup is very similar to the very recent work [9], making such an experiment feasible.