# Deep learning extraction of band structure parameters from density of
states: a case study on trilayer graphene

Paul Henderson  School of Computing Science, University of Glasgow, Scotland
Areg Ghazaryan  Institute of Science and Technology Austria (ISTA), Am Campus
1, 3400 Klosterneuburg, Austria  Alexander A. Zibrov  Department of Physics,
Harvard University, Cambridge, MA, USA  Department of Physics, University of
California, Santa Barbara, CA, USA  Andrea F. Young  Department of Physics,
University of California, Santa Barbara, CA, USA  Maksym Serbyn  Institute of
Science and Technology Austria (ISTA), Am Campus 1, 3400 Klosterneuburg,
Austria

###### Abstract

The development of two-dimensional materials has resulted in a diverse range
of novel, high-quality compounds with increasing complexity. A key requirement
for a comprehensive quantitative theory is the accurate determination of these
materials’ band structure parameters. However, this task is challenging due to
the intricate band structures and the indirect nature of experimental probes.
In this work, we introduce a general framework to derive band structure
parameters from experimental data using deep neural networks. We applied our
method to the penetration field capacitance measurement of trilayer graphene,
an effective probe of its density of states. First, we demonstrate that a
trained deep network gives accurate predictions for the penetration field
capacitance as a function of tight-binding parameters. Next, we use the fast
and accurate predictions from the trained network to automatically determine
tight-binding parameters directly from experimental data, with extracted
parameters being in a good agreement with values in the literature. We
conclude by discussing potential applications of our method to other materials
and experimental techniques beyond penetration field capacitance.

##  I Introduction

Electronic band structure of crystalline solids represents a simple yet very
rich example of emergence. Under the influence of scattering from the lattice
potential, the electron may acquire a different value of effective mass,
became massless, and acquire additional quantum numbers such as pseudospin. In
addition, electronic band structure determines basic properties of materials,
provided interactions are weak enough [1]. Therefore, identifying material
parameters that determine the band structure is of crucial importance. From a
theoretical point of view ab initio methods such as density functional theory
have achieved enormous success in this direction [2]. Nevertheless one
typically relies on experimental data to quantitatively extract band structure
properties. Experimentally there exist numerous ways to access the electronic
structure, such as angle resolved photoemission [3] and X-ray absorption
spectroscopy [4], de Haas-van Alphen effect based on magnetic oscillations
[5], analyzing reflection and absorption spectra [6], and electronic transport
measurement [7], to name just a few. Despite such a wealth of measurement
techniques, matching experimental results with theoretical predictions remains
a challenging problem due to the complexity of the band structure, which
translates into a large number of involved parameters.

The recent surge of two-dimensional (2D) materials [8] brings new aspects to
the problem of determining band structure. First, often the complexity and the
number of parameters in 2D materials is considerably lower compared to their
three-dimensional counterparts. Besides, 2D materials feature additional level
control such as modifying charge density by gating, and may have an extremely
high crystal quality. This opens access to high resolution experimental data,
which may potentially be used for precise determination of band structure
parameters. A particular example of such data is provided by so-called
penetration field capacitance measurements [9], that effectively probe the
density of states (DOS) of the material as a function of carrier density and
transverse electric field. Such experimental data has been used to determine
material parameters such as hopping matrix elements in several 2D systems [10,
11].

Typically, extraction of band structure parameters based on experimental data
relies on an efficient solution to what we term the _forward problem_. In the
specific example of penetration field capacitance measurements sensitive to
the DOS, this means simulating the DOS for specific values of material
parameters such as hopping matrix elements entering tight-binding model of the
band structure. However, the existence of an efficient solution for this
forward problem does not guarantee a fast solution to the inverse
problem—identification of the physical parameters corresponding to a set of
empirical data. The inverse problem is challenging because (i) solving for the
best-fitting parameters is a high-dimensional optimization problem that
requires numerous simulations of the forward problem at each step that can
quickly become very costly numerically; and (ii) experimental measurements are
typically affected by additional factors not easily accounted for in
simulation (e.g. geometric and parasitic capacitance, disorder), meaning that
an exact match between the data simulated in the forward problem and that
obtained from experiments is not possible. The typical approach is therefore
manual comparison of an experimental dataset with a large number of simulated
ones, relying on physical intuition of which features are important. This
process is laborious and computationally expensive 111This approach is aptly
referred to as graduate-student descent., calling for the development of more
efficient and systematic approaches.

In this work we present a machine learning based method that automates the
process of comparing numerical simulation and experimental data, so the
physical parameters of the band structure that gave rise to a particular
experimental dataset can be determined with minimal human effort. Recently
machine learning and artificial neural network techniques have seen various
applications in the realm of physical sciences [13]. In condensed matter
physics, artificial neural networks have been used to represent quantum states
[14, 15] and learn these states from available data [16, 17]. In a different
direction, recently machine learning models were use for photonic crystals
band diagram prediction and gap optimisation [18, 19, 20]. Despite a large
number of more theoretical applications, machine learning approaches are only
starting to be employed in analysis of experimental data. Recent examples
include identification of quantum phase transitions [21] and hidden orders
from experimental images [22]. These few examples highlight the strong
potential of machine learning based approaches on experimental data, that we
further exploit in the present work.

A conceptual overview of our approach is shown in Fig. 1. To extract the band
structure parameters from experimental data, we first train a deep neural
network (DNN) [23, 24] that solves the forward problem by replicating the
numerical calculation of the DOS (Section II.1). To this end we use the
simulation of the experimental data shown in Fig. 1(a). In the particular
example of penetration field capacitance data considered here, the simulator
uses the band structure parameters, the asymmetry potential between two edges
of the system (physically equivalent to transverse electric field) and the
chemical potential as input parameters. As an output we get charge density and
from that determine the DOS by differentiating density with respect to
chemical potential. A set of simulated data is used to train the DNN in Fig.
1(b). Constructed in a way to efficiently replace the data simulator, the DNN
acts as a function that takes the band structure parameters, the asymmetry
potential and directly charge density as input, and outputs the corresponding
DOS. It is constructed by learning from a large dataset of simulation results,
optimising its output to always match that of the simulator. The resulting DNN
represents a fast and differentiable replacement for the physical simulation.
It can therefore be used to efficiently solve the inverse problem (Section
II.2). In particular, the values of parameters that gave rise to a given
dataset are extracted using gradient-based optimisation in Fig. 1(c), where we
iteratively modify the band structure parameters until the DNN’s output
matches the provided DOS values.

Figure 1:  Our approach to the direct extraction of the band structure
parameters uses datasets obtained with a physical simulator shown in panel (a)
to train a DNN in panel (b), thus providing a more efficient solution to the
forward problem. The DNN-based simulator is used for a gradient optimization
of the band structure parameters, allowing to extract their values from
experimental data in panel (c).

The task of mapping a vector of inputs (e.g. band structure parameters) to a
continuous output (e.g. DOS) is known as regression in machine learning [24].
A DNN implements such a mapping as a series of chained matrix multiplications
(‘layers’) interleaved with elementwise non-linear functions (‘activations’).
Each layer multiplies the vector of outputs from the previous by some weight
matrix, to give an updated vector [23]; the final layer typical yields a
single value. The weight matrices are optimized (‘trained’) using first-order
optimization (e.g. gradient descent), such that the overall mapping from
inputs to output approximates some function defined by a training dataset of
inputs and desired outputs. The celebrated universal approximation theorem
[25] proves that a neural network with just two layers (but unbounded width)
can represent any smooth function. More recently, it has been shown that this
is also true for a neural network of bounded width (but unbounded depth) [26].
In practice, even finite-sized DNNs have proven very successful in
approximating complex functions in many domains of science. One of our
contributions is to show a DNN can also provide an accurate estimate of DOS
given band-structure parameters, field strength, and chemical potential.

Since our final goal is to determine the band-structure parameters from
experimental measurements of penetration capacitance, it might seem natural to
train the DNN for exactly this task (the inverse problem), instead of the
forward problem. However, this is not possible in practice. We have only a
_single_ experimental dataset, for which the parameters are unknown, whereas
machine learning techniques require a large dataset of training examples (with
the true output known) to learn the desired mapping from. If instead we
trained on easy-to-obtain simulated data, the resulting model would not work
on experimental data since the latter is significantly different from the
former both in terms of the relative magnitude of features in the data, and
the locations of features such as edges. These differences may arise since the
simulator uses a simplified effective model of the material and does not
account for screening at the microscopic level, disorder, strain, experimental
uncertainties, and possibly other ingredients. In contrast, statistical
learning theory requires that the training and test data be drawn from the
same distributions if a trained model is to work on the latter [27].

As a specific application, we demonstrate the framework outlined above on
Bernal stacked ($ABA$) trilayer graphene. For this material both band
structure parametrization [28, 29] and high quality experimental measurements
are readily available [11], calling for an accurate extraction of the band
structure parameters. The determination of the band structure was performed by
tour de force manual fitting in an earlier work [11], thus allowing us to
benchmark our approach. First, we train the DNN and show it gives an efficient
and accurate surrogate for numerical calculation of the DOS for this system,
for a wide range of band structure parameters (Section III.3). Next, we use
the DNN for automatically solving the inverse problem of determining the
physical parameters giving rise to certain values of the penetration
capacitance (Section III.4). Finally, we apply this to experimental data from
Ref. [11] by exploiting techniques from computer vision, that allow matching
important features of the measurements (e.g. Van Hove singularity peaks and
jumps of DOS), while ignoring features that differ between experimental and
simulated data (e.g. measurement noise and discreteness of calculation grid)
(Section III.5). The resulting values of parameters agree within error bars
with the estimates from the literature, thus providing a particular benchmark
for our approach.

The paper is organized as follows: Section II describes the general structure
of the DNN for the forward problem and our minimization approach for inverse
problem. Section III applies the framework constructed in Section II to $ABA$
graphene simulation results for the forward and inverse problems, eventually
utilizing it for extracting band parameters from experimental results. Finally
Section IV is devoted to discussion of the main results and generalization of
the model to other systems.

##  II Method

We assume access to a simulator $\mathcal{S}$ (see Fig. 1) that in the
particular example of the penetration field capacitance calculates the charge
density $n$ and density of states $\nu=\frac{\partial n}{\partial\mu}$ given
band structure parameters $\bm{\gamma}$, interlayer asymmetry $\Delta_{1}$,
and chemical potential $\mu$. We assume $\bm{\gamma}\in\Gamma$, where $\Gamma$
defines a physically-plausible range for those parameters, and similarly that
$(\Delta_{1},\,n)\in P$. The approach presented here is general, while the
specific physical meaning of these parameters will be discussed in Section
III. Typically the simulator will be slow to evaluate, making it difficult to
use ‘in the loop’ for solving the inverse problem, of finding the physical
parameters $\bm{\gamma}$ corresponding to observed data. We shall instead use
$\mathcal{S}$ to generate training data represented by tuples
$(\bm{\gamma},\,\Delta_{1},\,n,\,\nu)$ for a machine learning regression
model—a deep neural network—$F_{\omega}$, that will be trained to approximate
$\mathcal{S}$ in Section II.1. We shall then use the resulting DNN
$F_{\omega}$ when solving the inverse problem in Section II.2.

###  II.1 Forward problem

We introduce a function $F_{\omega}$ that maps $\bm{\gamma}$, $\Delta_{1}$ and
$n$ to the DOS. We choose $F_{\omega}$ to be a deep neural network (DNN) [23],
with weights $\omega$; these weights are free parameters that determine the
function it represents. Our goal is that $F_{\omega}$ matches $\mathcal{S}$ as
closely as possible for all relevant values of $\bm{\gamma}$, $\Delta_{1}$ and
$\mu$, i.e. if the simulator returns $n$ and $\nu$ for given
$(\bm{\gamma},\,\Delta_{1},\,\mu)$ and if $(\Delta_{1},n)\in P$ are within the
domain of physically realistic parameters, then DNN approximates well the DOS,
$F_{\omega}(\bm{\gamma},\Delta_{1},n)\approx\nu$. The network weights $\omega$
would ideally be set to minimise the absolute difference between the network’s
predicted values and those $\nu$ from the simulator, over the entire parameter
space $\Gamma\times P$:

$\omega=\mathop{\rm argmin}{\omega^{\prime}}\int_{(\Delta_{1},\,n)\in
P,\bm{\gamma}\in\Gamma}\big{|}F_{\omega^{\prime}}(\bm{\gamma},\,\Delta_{1},\,n)-\\\
\nu\big{|}\,\mathrm{d}\bm{\gamma}\,\mathrm{d}\Delta_{1}\mathrm{d}n$ (1)

In practice we instead minimise the mean error over a finite training set [24]
of points $\mathcal{T}\subset\Gamma\times P$ at which we have precomputed $n$
and $\nu$ using the simulator $\mathcal{S}$, i.e.

$\omega=\mathop{\rm
argmin}{\omega^{\prime}}\;\frac{1}{|\mathcal{T}|}\sum_{(\bm{\gamma},\,\Delta_{1},\,n,\,\nu)\in\mathcal{T}}\big{|}F_{\omega^{\prime}}(\bm{\gamma},\,\Delta_{1},\,n)-\nu\big{|}.$
(2)

To solve this optimisation problem, the weights $\omega$ are initialised using
the heuristic of Ref. [30], then iteratively updated using the first-order
stochastic-gradient optimiser Adam [31] with a minibatch size of 512 and
learning rate (step size) of $10^{-3}$. We use a DNN with five fully-connected
layers of 512 units each, with ELU nonlinearities [32], layer normalisation
[33], and residual connections [34]. For the input layer, we use Fourier
feature embedding with four octaves [35]; for the output layer, we use a
single linear unit, see Appendix A for details. We select these architectural
parameters, and determine when to stop training the DNN, based on its
performance on a separate validation set, that is disjoint from $\mathcal{T}$.
We implemented the DNN using the TensorFlow library [36]. Codes for generating
data, training, and evaluating the DNN are publicly available [37].

###  II.2 Inverse problem

Suppose we have an experimental dataset that provides measurements of
penetration capacitance $C_{p}(D,\,n)$, which is a quantity that is sensitive
to the DOS, as we discuss below. The measurements are acquired while varying
$(D,\,n)$ over some finite set $Q$; here $D$ is the strength of an externally-
applied electric field that affects $\Delta_{1}$. The inverse problem is to
determine the band structure parameters $\bm{\gamma}$ for the system, i.e. to
find the setting for $\bm{\gamma}$ for which the simulated results are the
closest to the experimental ones. We use the trained DNN $F_{\omega}$ in place
of the simulator $\mathcal{S}$; we therefore seek

$\bm{\gamma}=\mathop{\rm
argmin}_{\bm{\gamma}^{\prime}\in\Gamma}\,\min_{\bm{\alpha},\,\bm{\beta}}\,\frac{1}{|Q|}\sum_{(D,n)\in
Q}d\\!\left(C_{p}^{*}(D,n),\,C_{p}(D,n)\right)$ (3)

where $d$ is a metric quantifying the difference between experimental and
predicted values, and we introduced a function $C_{p}^{*}(D,n)$ that relies on
the electrostatic model to convert the DOS approximated by the DNN into the
penetration field capacitance value,

$C_{p}^{*}(D,n)=\frac{\beta_{1}}{\beta_{2}+F_{\omega}(\bm{\gamma}^{\prime},\,\alpha_{1}D+\alpha_{3}D^{3},n)},$
(4)

with $\alpha_{1,3}$ and $\beta_{1,2}$ being parameters that encode the
screening and electrostatic characteristics of the experimental setup
respectively. This function and intuition behind these additional parameters
is described in detail in Section III.1.

We solve the minimization problem (3) using the Adam optimiser [31] 222We also
experimented with the constrained trust-region method of Ref. [52], and the
quasi-Newton method L-BFGS-B [53], both as implemented in scipy [54]. However,
noise in the experimental data caused both these methods to perform poorly,
becoming trapped in local optima.. This optimiser makes use of the gradient of
the objective with respect to $\bm{\gamma}$; to calculate this, we use
reverse-mode automatic differentiation on the objective in Eq. (3), similar to
the back-propagation process used when training the DNN [39], but now with the
weights $\omega$ held fixed and gradients instead propagated to the inputs
$\bm{\gamma}$.

##  III Application to ABA graphene

In this section we discuss a specific application of our method—determining
the band structure of Bernal stacked trilayer (“$ABA$”) graphene. We first
describe the physical system and how it is simulated in Section III.1. Then,
we discuss the dataset that we generate using this simulator (Section III.2)
and use to train the DNN. Finally, we discuss the performance of the trained
DNN when solving the forward (Section III.3) and inverse (Section III.4)
problems, and how we apply it to experimental data (Section III.5).

###  III.1 Physical system and simulation

The band structure of $ABA$ graphene can be decomposed through rotation of the
basis into monolayer and bilayer graphene type sectors, which get coupled
through a displacement field applied between the layers. The Hamiltonian
matrix takes the form [29]

$H_{\mathbf{k}}=\left(\begin{array}[]{cc}H_{\mathrm{SLG}}&V\\\
V^{T}&H_{\mathrm{BLG}}\end{array}\right),$ (5)

where

$\displaystyle H_{\mathrm{SLG}}$
$\displaystyle=\left(\begin{array}[]{cc}\Delta_{2}-\frac{\gamma_{2}}{2}&v_{0}\pi^{\dagger}\\\
v_{0}\pi&-\frac{\gamma_{5}}{2}+\delta+\Delta_{2}\end{array}\right),$ (8)
$\displaystyle H_{\mathrm{BLG}}$
$\displaystyle=\left(\begin{array}[]{cccc}\frac{\gamma_{2}}{2}+\Delta_{2}&\sqrt{2}v_{3}\pi&-\sqrt{2}v_{4}\pi^{\dagger}&v_{0}\pi^{\dagger}\\\
\sqrt{2}v_{3}\pi^{\dagger}&-2\Delta_{2}&v_{0}\pi&-\sqrt{2}v_{4}\pi\\\
-\sqrt{2}v_{4}\pi&v_{0}\pi^{\dagger}&\delta-2\Delta_{2}&\sqrt{2}\gamma_{1}\\\
v_{0}\pi&-\sqrt{2}v_{4}\pi^{\dagger}&\sqrt{2}\gamma_{1}&\frac{\gamma_{5}}{2}+\delta+\Delta_{2}\end{array}\right),$
(13) $\displaystyle V$
$\displaystyle=\left(\begin{array}[]{cccc}\Delta_{1}&0&0&0\\\
0&0&0&\Delta_{1}\end{array}\right).$ (16)

Here $v_{i}=\sqrt{3}a\gamma_{i}/2$ and $a=2.46\,\AA$ is the lattice constant.
$\gamma_{i}$ and $\delta$ are the hopping and on-site potential parameters of
the physical system. $\pi=\xi k_{x}+ik_{y}$ is the crystal imaginary momentum
measured with respect $K^{\pm}$ valley points labeled by $\xi=\pm 1$.
$\Delta_{1}$ describes the effect of the external field and $\Delta_{2}$
charge asymmetry between external and internal layers.

For producing simulated data we calculate density and DOS on a grid using
eigenvalues of Hamiltonian (5) $\epsilon_{\mathbf{k}}$ for each value of
chemical potential $\mu$

$\displaystyle n\left(\mu\right)$
$\displaystyle=4\frac{S_{k}}{\left(2\pi\right)^{2}}\frac{1}{N}\sum_{\mathbf{k}}n_{F}\left(\epsilon_{\mathbf{k}}-\mu\right),$
(17) $\displaystyle\nu\left(\mu\right)$
$\displaystyle=4\frac{S_{k}}{\left(2\pi\right)^{2}}\frac{1}{N}\sum_{\mathbf{k}}n^{\prime}_{F}\left(\epsilon_{\mathbf{k}}-\mu\right),$
(18)

where $4$ accounts for spin and valley degeneracy, $S_{k}$ is the area of the
Brillouin zone sampled in the grid, $N$ is the number of grid points and
$n_{F}\left(x\right)=1/\left(e^{x/T}+1\right)$ is the Fermi-Dirac
distribution. We use finite temperature $T=0.025\,\mathrm{meV}$ in the
calculation to smoothen singularities of DOS at Van Hove singularity points,
and a hexagonal grid with a momentum cutoff $ka=0.15$ and $7.5\times 10^{5}$
grid points.

To obtain the dependence of the experimentally measured penetration field
capacitance on the DOS we can use the following formula [40, 41, 11]:

$C_{p}=\frac{c_{t}c_{b}}{c_{t}+c_{b}+e^{2}\nu},$ (19)

where $c_{b}$ and $c_{t}$ is the geometric capacitance of bottom and top gate,
$e$ is the charge of the electron. However, this is a simplified formula which
ignores layer polarization change due to the applied electric field $D$, and
views the trilayer graphene as a single layer system. Therefore, to allow a
more general relationship, we introduced the parameters $\beta_{1}$ and
$\beta_{2}$ in (4). Besides linear screening, we allow for non-linear
component, so

$\Delta_{1}=\alpha_{1}D+\alpha_{3}D^{3}.$ (20)

Due to the presence of a parasitic capacitance, we model the relation between
experimentally measured and calculated charge density by

$n^{\prime}=n+\left(\mu-\mu\left(\Delta_{1}=0,n=0\right)\right)/\eta,$ (21)

where $\eta$ is the inverse of parasitic capacitance.

In the figures, simulation results are presented in terms of inverse of DOS
$\nu^{-1}$ in units of $\mathrm{eV}\cdot A_{\mathrm{u.c}}$, where
$A_{\mathrm{u.c}}=\sqrt{3}a^{2}/2$ is the area of graphene unit cell. Charge
density is always presented in units of $\mathrm{cm^{-2}}$. Band structure
parameters $\gamma_{i}$ and parameters $\delta$, $\Delta_{1}$, $\Delta_{2}$
and $\mu$ are presented in $\mathrm{meV}$.

Figure 2:  Predictions of the inverse DOS, $\nu^{-1}$, by our DNN, compared
with results from a numerical simulator. (a) and (b) corresponds to the band-
structure parameters $\gamma_{2}=-4.3$, $\gamma_{3}=-342$, $\gamma_{4}=133.4$,
$\Delta_{2}=-0.4$ and $\gamma_{2}=-15.6$, $\gamma_{3}=-248.8$,
$\gamma_{4}=104.9$, $\Delta_{2}=11.2$ respectively. Both of the sets were not
present in the training data for the DNN. For each we show: (left) $\nu^{-1}$
as a function of field strength $\Delta_{1}$ and electron charge density $n$,
as predicted by a simulator; (middle) $\nu^{-1}$ as predicted by our DNN;
ideally this would exactly match the simulator output; (right) the absolute
difference between the simulator and DNN values. We see good correspondence
between the DNN’s outputs and the simulator, indicating that it has been
successfully trained to replicate the simulator output, at much lower
computational cost.

###  III.2 Data generation

Table 1: Range of physical parameter and grid step used in the calculations. Parameter |  Minimum (meV) |  Maximum (meV) |  Step (meV)
---|---|---|---
$\gamma_{0}$ |  $3100$ |  $3100$ |  –
$\gamma_{1}$ |  $380$ |  $380$ |  –
$\gamma_{2}$ |  $-25$ |  $-1$ |  $3$
$\gamma_{3}$ |  $-340$ |  $-242$ |  $7$
$\gamma_{4}$ |  $100$ |  $177$ |  $7$
$\gamma_{5}$ |  50 |  50 |  –
$\delta$ |  35.5 |  35.5 |  –
$\Delta_{2}$ |  $0$ |  $11$ |  $1$

Using the simulation $\mathcal{S}$ described above, we generate data to train
and evaluate the DNN. In order to define a suitable grid, we choose the space
$\Gamma$ of valid physical parameters using physical insights into the meaning
of the tight-binding parameters. Since $\delta$ and $\gamma_{5}$ introduce
only a simple energy shift on the spectrum for a considerable range of values,
we fix them to $\delta=35.5\,\mathrm{meV}$ and $\gamma_{5}=50\,\mathrm{meV}$.
We also fix $\gamma_{0}=3100\,\mathrm{meV}$ and
$\gamma_{1}=380\,\mathrm{meV}$; these are intralayer and direct interlayer
hopping matrix elements, which were studied extensively in graphite, and in
monolayer and bilayer graphene. The remaining tuning parameters that
constitute $\bm{\gamma}$ are $\gamma_{2}$, $\gamma_{3}$, $\gamma_{4}$,
$\Delta_{2}$ and $\eta$. These parameters are allowed to take values in the
ranges shown in Table 1. Finally, we define $P$, the region of acceptable
parameters $\Delta_{1}$ and $n$, according to $0\leq\Delta_{1}\leq
75\,\mathrm{meV}$ and $|n|<4\times 10^{12}\,\mathrm{cm^{-2}}$, corresponding
to the physically achievable range of transverse electric fields and densities
respectively.

To generate datasets for training, validation, and testing of the DNN, we
iterate over a grid of $\gamma_{2}$, $\gamma_{3}$, $\gamma_{4}$ and
$\Delta_{2}$, with the respective increments for each specified in the last
column of Table 1. The increment is chosen so that for each parameter the grid
contains on the order of 10 points. We also add an additional random offset of
up to half the step size to each value, to avoid the DNN overfitting on points
of the grid. We then partition this set of parameter vectors into three
disjoint subsets: we randomly assign 75% of parameter values to the training
set $\mathcal{T}$, 10% for validation and the remainder for testing 333Note
that the validation set is used for determining convergence and selecting the
network architecture, while all results are reported on the testing set; this
avoids overfitting to the testing data..

For each value of the structure parameters, we vary $\Delta_{1}$ and $\mu$ in
the intervals $\left[0,75\right]\,\mathrm{meV}$ and
$\left[-115,115\right]\,\mathrm{meV}$ with steps $2.5\,\mathrm{meV}$ and
$0.025\,\mathrm{meV}$ respectively. We calculate the corresponding values of
$n$ and $\nu$ using $\mathcal{S}$ as specified in Eqs. (17)-(18). Since the
charge density $n$ for different values of the parasitic capacitance may be
determined straightforwardly from that when it is zero, we fix the inverse
parasitic capacitance $\eta=\infty$. This yields a dataset of tuples of values
$(\gamma_{2},\,\gamma_{3},\,\gamma_{4},\,\Delta_{2},\,\Delta_{1},\,n,\,\nu)$.
We discard tuples where $n\notin P$ (which is beyond our range of interest),
or both $|n|<0.01\times 10^{12}\,\mathrm{cm^{-2}}$ and
$\nu<0.01\mathrm{eV}\cdot A_{\mathrm{u.c}}$ (i.e. very close to the charge
neutrality point and the gapped region, respectively). During the training
process itself, we sample the value of inverse parasitic capacitance, $\eta$,
from a log-uniform distribution on the interval
$\left[50,5000\right]\,\mathrm{10^{-12}\,meV\times cm^{2}}$, independently for
each element of each minibatch, and transform $n$ according to (21). In total
our training set contains $1.38\times 10^{9}$ points from 14580 settings of
the band structure parameters; there are $184\times 10^{6}$ points in the
validation set, and $276\times 10^{6}$ in the test set.

###  III.3 Forward problem: calculating $\nu$

Using the dataset of simulation results for $ABA$ graphene described above, we
train a DNN following the protocol in Section II.1, to predict the density of
states $\nu$.

To evaluate the accuracy of the trained DNN, we use it to predict $\nu$ for
points in the held-out test set, i.e. for parameters
$(\gamma_{2},\gamma_{3},\gamma_{4},\Delta_{2})\notin\mathcal{T}$ that were not
contained in the training dataset. Quantitatively, our method achieves a mean
absolute error in $\nu$ of $0.00288\,\mathrm{eV}\cdot A_{\mathrm{u.c}}$,
corresponding to a mean relative error of $5.95\%$. We further analyzed the
cause of these errors, by resimulating a random subset 300 of plots in the
test set using a finer momentum grid with $6.75\times 10^{6}$ points and
cutoff $ka=0.2$. Comparing the DNN’s predictions to these more accurate
simulations, the absolute error reduces to $0.00067\,\mathrm{eV}\cdot
A_{\mathrm{u.c}}$, corresponding to a mean relative error of $1.32\%$. We also
measure the relative error of the lower-resolution simulations forming our
dataset, with respect to the 300 higher-resolution plots, and find this is
$5.97\%$. Remarkably, the DNN’s predictions are therefore slightly more
accurate than the lower-resolution simulations used for its training data. We
hypothesize that this is due to the DNN’s inductive bias towards learning
smooth functions discouraging it from learning the very variable high-
frequency artifacts that arise due to the coarser momentum grid. Therefore, in
the following figures we show simulator results for higher resolution momentum
grid, although DNN training was performed using lower-resolution data.

Examples of outputs from the DNN, alongside corresponding outputs from the
simulator, are shown in Figure 2. We also visualize the mean absolute
difference between the DNN’s predictions, and those of the simulator. We see
that the DNN produces outputs closely matching the simulator, most importantly
accurately incorporating features such as discontinuities and (smoothened)
singularities that originate from the changes of Fermi surface topology upon
tuning density, $n$ or $\Delta_{1}$. For instance, the diagonal feature for
negative densities in Fig. 2(a) corresponds to the disappearance of the hole
pocket upon increasing density, while parabola-shaped blue line in the
negative density region corresponds to the Van Hove singularity where DOS
would show a logarithmic divergence in absence of cutoffs introduced due to
finite temperature and finite grid size in momentum space. The plots of
absolute differences show that the errors are typically very small in regions
of uniform $\nu^{-1}$, and are dominated by small misalignments of the
features that correspond to DOS discontinuities or singularities.

We also measure the difference in computation time for the simulator and the
DNN. Running with four threads on a 3.6GHz Intel Core i7-9700K CPU, the
simulator takes 1445 s to evaluate $\nu$ for the grid of $(\Delta_{1},\mu)$
specified in Section III.2 and the fine momentum grid. Running the DNN on the
same CPU with four threads takes just 1.1 s. Moreover, it does so with higher-
resolution sampling of the range of values of $n$ that are of interest, due to
taking $n$ as an input instead of $\mu$, hence avoiding the need to discard
calculations where $n$ is outside the range of interest.

| $\gamma_{2}$ | $\gamma_{3}$ | $\gamma_{4}$ | $\Delta_{2}$
---|---|---|---|---
true | $-15.58$ | $-314.1$ | 161.1 | $-0.15$
predicted | $-15.36$ | $-313.3$ | 160.2 | $-0.16$

| $\gamma_{2}$ | $\gamma_{3}$ | $\gamma_{4}$ | $\Delta_{2}$
---|---|---|---|---
true | $-18.64$ | $-265.7$ | 114.8 | 5.30
predicted | $-18.56$ | $-265.0$ | 115.1 | 5.28

Figure 3: Results on inverting simulator outputs, i.e. predicting the band-structure parameters that gave rise to a given plot. Each panel shows (left) simulated $\nu^{-1}$ as a function of $\Delta_{1}$ and $n$, for the specified ‘true’ band parameters; (center) simulated $\nu^{-1}$ for the values of the band parameters inferred by our method; (right) absolute difference between simulations and inferred data. The tables at the top compare the true and predicted values of parameters.  Table 2: Quantitative performance on recovering band-structure parameters from simulated data. Parameter |  Mean absolute error (meV) |  Median abs. error (meV) |  Median rel. error
---|---|---|---
$\gamma_{2}$ |  $0.12$ |  $0.10$ |  $0.0156$
$\gamma_{3}$ |  $1.54$ |  $1.30$ |  $0.0512$
$\gamma_{4}$ |  $0.59$ |  $0.45$ |  $0.0240$
$\Delta_{2}$ |  $0.04$ |  $0.03$ |  $0.0135$

###  III.4 Inverse problem: determining parameters

We next consider the inverse problem, of finding the band-structure parameters
corresponding to a set of measurements. We first evaluate our method on
simulated data; this allows validating the proposed approach quantitatively,
by measuring how accurately it recovers the parameters that were input to the
simulator. Specifically, we select 100 simulator-generated plots at random
from the held-out test set to use as input; each plot shows how DOS $\nu$
varies with $\Delta_{1}$ and $n$, and we aim to find the corresponding band-
structure parameters.

To find the parameters $\bm{\gamma}$, we follow the procedure described in
Section II.2, using the same DNN as discussed in Section III.3. As the inputs
are simulated data, the DOS itself is directly available, hence we choose the
metric $d$ in (3) to be the mean squared difference between the input and
predicted $\nu$.

Table 2 shows quantitative results from our approach, calculated over 100
plots from the test set. For each band-structure parameter, we give the mean
and median absolute error, and also the median relative error after first
subtracting the midpoint of the corresponding range in Table 1. We see that
both absolute values of error and also its relative value are small, with the
greatest relative error being $5\%$ for the parameter $\gamma_{3}$, and even
smaller for the remaining parameters. Such a confidence interval in
determining tight binding parameters is much smaller compared to the error
bars typically available in the literature. For instance, Ref. [11] determines
$\gamma_{3}$ parameter with the precision of $17\%$. Further discussion on
quantitative performance of the model using parity plots is provided in
Appendix B. Qualitative illustrations of our results for two particular points
in the parameter space are given in Figure 3. We show the plot provided as
input to our method, and the result of running the simulator on the predicted
parameter values. Here we show high-resolution momentum grid images for both
the input and the prediction, although the DNN operated on low-resolution
momentum grid images. We see that the residual error is very small, and is
less concentrated in the vicinity of DOS discontinuities or singularities,
compared to Fig. 2 that benchmarked the forward approach.

Figure 4: Inferring band-structure parameters from experimental data. The left
plot shows experimental $C_{p}$ measurement as a function of displacement
field $D$ and charge density $n$; the middle and right plots show the inverse
of DOS $\nu^{-1}$ as a function of $\Delta_{1}$ and $n$ for parameters
$\gamma_{2}=-21$, $\gamma_{3}=-290$, $\gamma_{4}=141$, $\Delta_{2}=3.5$ and
$\gamma_{2}=-15.5$, $\gamma_{3}=-312.9$, $\gamma_{4}=132.3$, $\Delta_{2}=3.0$
respectively. The middle figure is the best match to the experimental results
cited in the literature [11], while the right one is the prediction of the
DNN. The middle and right plots were obtained using simulator with higher
resolution of $\mu$ ($\Delta\mu=0.002$), hexagonal grid with maximum momentum
$ka=0.2$, $6.75\times 10^{6}$ grid points, larger range of $\Delta_{1}$ and
with a step of $0.5$.

###  III.5 Application to experimental data

Finally, we apply our method to the high quality experimental dataset of
penetration capacitance measurements for $ABA$ graphene used in Ref. [11].
Note, that Ref. [11] utilized additional datasets, namely, the density of
states in presence of magnetic field (Landau level regime) in order to assist
in determination of tight binding parameters and to infer the associated error
bars. In contrast, here we use solely zero magnetic field penetration field
capacitance data. Our goal is to solve a similar inverse problem to the one in
Section III.4, i.e. predicting the band-structure parameters for this system,
but now given only noisy measurements of $C_{p}$.

Using the mean squared difference between predicted and experimental $C_{p}$
for the metric $d$ in Eq. (3) yields poor results in this setting. This is
because features in the simulated and experimental plots (e.g. jumps in DOS)
do not perfectly align for any choice of parameters, and also because there is
a small residual difference in the absolute values of $C_{p}$ between
simulated and experimental data even when the plots are optimally aligned.
These discrepancies can be at least partially attributed to multiple physical
mechanisms that are not incorporated in our model. In particular, we take
account of the screening only phenomenologically by using the third order
expansion (20), whereas the use of self-consistent Hartree-Fock screening
would result in a more complex function. In addition, we ignore potential
renormalization of the band structure parameters by the strong applied
electric fields, that could make all band structure parameters dependent on
$n$ and $\Delta_{1}$. Finally, disorder and strain effects are not
incorporated in our model, yet these may lead to inhomogeneous broadening of
the features.

To mitigate this issue, we adopt techniques from computer vision that are used
to align dissimilar images; these are typically employed in tasks such as
template matching [43, 44], medical image registration [45, 46], and alignment
of satellite imagery [47]. Firstly, in order to match prominent features in
the plot such as edges regardless of the absolute values of $C_{p}$, we match
the derivative of $C_{p}$ with respect to $D$ and $n$ instead of $C_{p}$
itself. Secondly, instead of the mean squared error, we use the negative zero-
normalised cross-correlation [44], which is invariant to changes in mean and
standard deviation. We therefore set

$d(C^{*}_{p},\,C_{p})=-\frac{\partial_{D}C^{*}_{p}-m^{*}_{D}}{\sqrt{\left\langle(\partial_{D}C^{*}_{p}-m^{*}_{D})^{2}\right\rangle}}\cdot\frac{\partial_{D}C_{p}-m_{D}}{\sqrt{\left\langle(\partial_{D}C_{p}-m_{D})^{2}\right\rangle}}\\\
-\frac{\partial_{n}C^{*}_{p}-m^{*}_{n}}{\sqrt{\left\langle(\partial_{n}C^{*}_{p}-m^{*}_{n})^{2}\right\rangle}}\cdot\frac{\partial_{n}C_{p}-m_{n}}{\sqrt{\left\langle(\partial_{n}C_{p}-m_{n})^{2}\right\rangle}}$
(22)

in Eq. (3), where $\left\langle\,\cdot\,\right\rangle$ is the mean over all
$(D,\,n)$, $m^{*}_{D}=\left\langle\partial_{D}C^{*}_{p}\right\rangle$,
$m^{*}_{n}=\left\langle\partial_{n}C^{*}_{p}\right\rangle$,
$m_{D}=\left\langle\partial_{D}C_{p}\right\rangle$,
$m_{n}=\left\langle\partial_{n}C_{p}\right\rangle$, and the derivatives of the
experimental measurements are approximated using finite differences. Lastly,
in order to mitigate against local minima, we average the objective (22) over
multiple scales, i.e. subsampling the grid $Q$ by factors $2^{i}$ for
$i\in\\{1\ldots 5\\}$, and correspondingly smoothing and downsampling $C_{p}$
and $C_{p}^{*}$ to form scale-space pyramids [48].

Solving the resulting minimization problem as described in Section II.2
recovers the following parameters:

$\displaystyle\gamma_{2}$ $\displaystyle=$
$\displaystyle-15.5\,\mathrm{meV},\quad\gamma_{3}=-312.9\,\mathrm{meV},$
$\displaystyle\gamma_{4}$ $\displaystyle=$ $\displaystyle
132.3\,\mathrm{meV},\quad\Delta_{2}=3.0\,\mathrm{meV}.$ (23)

These are broadly consistent with the estimates obtained from manual fitting
in Ref. [11]. More precisely, our estimates of $\gamma_{2}$, $\gamma_{3}$ and
$\gamma_{4}$ lie within the error bars of [11] (although $\gamma_{2}$ is at
the upper limit), while $\Delta_{2}$ is slightly outside (we predict 3.0,
compared with their confidence interval $[3.3,\,3.7]$). Notably, the obtained
$\gamma_{2}$ is actually very close to the value obtained for ABC trilayer
graphene [49].

In Figure 4, we visualize the experimental data, and the output from the
simulator when run with the above parameters at high resolution. The simulated
plot shows a good match to the experimental one, with similar features (e.g.
discontinuities) appearing at similar locations. We also show the result using
the parameters of Ref. [11], these are again qualitatively very similar. From
visual inspection the separation between the two Van Hove singularities in the
hole region appears to be wider for predicted parameters than for the one from
Ref. [11]. This trend is consistent with the experimental results. For all
other features the difference between the two models is visually
indistinguishable. Therefore, the newly obtained parameters are a better fit
to the experiment data.

Table 3:  Quantitative results of sensitivity analysis on experimental data. For each parameter, we record the values to which it must be decreased/increased to cause the distance function $d(C_{p}^{*},\,C_{p})$ to increase by $1\%$. We specify both absolute and relative values; a large percentage indicates low sensitivity to that parameter/direction, as a large change is required to affect the quality of fit.  Parameter | Decreasing | Best fit (meV) | Increasing
---|---|---|---
| meV | rel.diff. | meV | rel.diff.
$\gamma_{2}$ | $-17.25$ | $-11.2\%$ | $-15.51$ | $-14.45$ | $+6.8\%$
$\gamma_{3}$ | $-343.5$ | $-9.8\%$ | $-312.9$ | $-296.8$ | $+5.1\%$
$\gamma_{4}$ | $111.1$ | $-16.1\%$ | $132.3$ | $153.2$ | $+15.8\%$
$\Delta_{2}$ | $2.72$ | $-9.5\%$ | $3.00$ | $3.17$ | $+5.5\%$

We also measure the sensitivity of the experimental fit with respect to each
of the band structure parameters. Specifically, we increase and decrease each
of the parameters $\gamma_{2}$, $\gamma_{3}$, $\gamma_{4}$ and $\Delta_{2}$
separately, until the distance $d(C_{p}^{*},\,C_{p})$ between predicted and
experimental $C_{p}$ values increases by more than $1\%$ from the ‘best fit’
value. This lets us measure quantitatively how large a change in each
parameter is required to cause the same reduction in the quality of fit.
Results from this analysis are given in Table 3. We see that the fit is least
sensitive to $\gamma_{4}$, with changes of around $16\%$ in either direction
required to cause a change of $1\%$ in the distance $d$. Conversely,
increasing $\gamma_{3}$ by just $5.8\%$ leads to a $1\%$ change in $d$;
$\Delta_{2}$ and $\gamma_{2}$ are similarly sensitive to increases.

##  IV Discussion

To summarize, we proposed a deep neural network based approach to determining
band structure parameters of two-dimensional materials. Our framework consists
of two steps, that rely on the existence of a method to simulate the
experimentally accessible data. In the first step we train a deep neural
network to obtain a more efficient replacement for the data simulator. In the
second step, we extract band structure parameters by minimizing the difference
between the experimental dataset and the data simulated by the neural network.
To illustrate the application of our framework, we focused on a specific
material — trilayer graphene — and performed band structure parameter
extraction using the experimental data on the penetration field capacitance
that effectively probes the density of states. Our procedure resulted in the
precise determination of band structure parameters that are close to those
determined manually.

In contrast to manual fitting, the approach proposed in our work is more
automated, thus involving less human effort. In addition, it is capable of
providing more precise values of band structure parameters, and of estimating
associated error bars. We expect that our framework can be easily applied with
minimal modifications to the penetration field capacitance data in other two
dimensional materials. Moreover, it can be generalized to other experimental
probes, such as penetration field capacitance in presence of magnetic field,
transport, scanning tunneling microscopy, and other probes that are sensitive
to the band structure details. Application to different experimental
quantities requires replacement of the simulator function in our framework,
which would be straightforward. A more intricate step consists in
understanding the relation between the simulated quantities and experimental
measurements. For instance, in the case of penetration field capacitance, the
matching of experimental data additionally required incorporating non-trivial
screening of electric field, the presence of parasitic capacitance and
geometric capacitances of the gates. These steps require physical insight into
the details of the experimental measurements.

Finally, a conceptually novel generalization direction may include
incorporating interaction effects into our framework. For instance, the recent
observation of Stoner transitions in twisted bilayer graphene [50, 51] or
chirally stacked $ABC$ trilayer graphene [49], calls for unambiguous
identification of interaction parameters alongside the band structure details
(that may be in turn quite significantly renormalized due to interactions). We
expect that such an extension of our framework may be feasible, provided one
is able to construct an efficient simulator of the density of states that
incorporates the interactions. This would constitute an important step towards
an experimental extraction of the complete Hamiltonian governing electronic
degrees of freedom, thus bringing two dimensional materials on par with
quantum simulators that use cold atoms or other platforms to synthetically
engineer and verify a desired Hamiltonian.

## Acknowledgments

A.F.Y. acknowledges primary support from the Department of Energy under award
DE-SC0020043, and additional support from the Gordon and Betty Moore
Foundation under award GBMF9471 for group operations.

\begin{overpic}[width=433.62pt,abs,unit=1mm]{dnn-arch}
\put(1.0,21.0){$\gamma_{2,3,4},\Delta_{2}$}
\put(3.0,17.0){$\eta,\Delta_{1},n$} \put(176.5,18.75){$\nu$} \end{overpic}
Figure 5: Architecture of the deep neural network $F_{\omega}$. The DNN takes
band structure parameters
$(\gamma_{2},\,\gamma_{3},\,\gamma_{4},\,\Delta_{2})$, inverse parasitic
capacitance $\eta$, interlayer asymmetry $\Delta_{1}$ and charge density $n$
as input. It processes them with a Fourier embedding layer (FE), linear layer
with 512 outputs, and ELU nonlinearity, followed by five blocks with residual
connections consisting of linear, ELU, and layer-norm (LN) layers. The
resulting features are mapped via a linear projection to the density of states
$\nu$.

##  Appendix A Additional model and training details

The architecture of our deep neural network $F_{\omega}$ is shown in Fig. 5.
In total, the DNN has $1.35\times 10^{6}$ weights (i.e. trainable parameters).
We train the model using Adam [31] with a fixed learning rate of $10^{-3}$,
and the standard momentum parameters $\beta_{1}=0.9$, $\beta_{2}=0.999$ and
$\epsilon=10^{-7}$. The batch size is set to 512. Convergence is determined
according to the loss on the held-out validation set – we stop training once
the validation loss has failed to decrease for 10 consecutive epochs, and
retain the checkpoint with minimum validation loss. We use a single Nvidia RTX
2080 Ti GPU for training, operating at 32-bit precision.

Figure 6: Parity plots of the band structure parameters from simulated data
for the inverse problem using the test set of 100 plots. The analyzed data is
the same as in Table 2. MAE (mean absolute error), RMSE (root mean square
error) and $R^{2}$ metrics are also shown. Equality line ($y=x$) is shown in
red.

##  Appendix B Additional quantitative results on the inverse problem

To gain further insights on the performance of the model for the inverse
problem, Fig. 6 shows parity plots of the band structure parameters
$\gamma_{2}$, $\gamma_{3}$, $\gamma_{4}$ and $\Delta_{2}$ for the same test
set of 100 plots which was used for obtaining the metrics in Table 2. As noted
in the discussion of Table 2 the largest deviation is observed for parameter
$\gamma_{3}$, and Fig. 6 (b) shows that the prediction of the model is always
higher than the true values. This systematic error hints at a direction for
further improving the model accuracy, and would be an interesting target for
future investigation.

## References

  * Harrison [1989] W. A. Harrison,  _Electronic structure and the properties of solids: the physics of the chemical bond_ (Dover Publications, New York, 1989).
  * Hasnip  _et al._ [2014] P. J. Hasnip, K. Refson, M. I. Probert, J. R. Yates, S. J. Clark, and C. J. Pickard, Density functional theory in the solid state, Phil. Trans. R. Soc. A 372, 20130270 (2014).
  * Damascelli  _et al._ [2003] A. Damascelli, Z. Hussain, and Z.-X. Shen, Angle-resolved photoemission studies of the cuprate superconductors, Rev. Mod. Phys. 75, 473 (2003).
  * De Groot [2001] F. De Groot, High-resolution x-ray emission and x-ray absorption spectroscopy, Chem. Rev. 101, 1779 (2001).
  * Shoenberg [2009] D. Shoenberg,  _Magnetic oscillations in metals_ (Cambridge university press, 2009).
  * Wooten [2013] F. Wooten,  _Optical properties of solids_ (Academic press, 2013).
  * Mizutani [2001] U. Mizutani,  _Introduction to the electron theory of metals_ (Cambridge university press, 2001).
  * Novoselov  _et al._ [2016] K. Novoselov, A. Mishchenko, A. Carvalho, and A. Castro Neto, 2D materials and van der Waals heterostructures, Science 353 (2016).
  * Eisenstein  _et al._ [1992] J. Eisenstein, L. Pfeiffer, and K. West, Negative compressibility of interacting two-dimensional electron and quasiparticle gases, Phys. Rev. Lett. 68, 674 (1992).
  * Island  _et al._ [2019] J. Island, X. Cui, C. Lewandowski, J. Khoo, E. Spanton, H. Zhou, D. Rhodes, J. Hone, T. Taniguchi, K. Watanabe, L. S. Levitov, M. P. Zaletel, and A. F. Young, Spin–orbit-driven band inversion in bilayer graphene by the van der Waals proximity effect, Nature 571, 85 (2019).
  * Zibrov  _et al._ [2018] A. A. Zibrov, P. Rao, C. Kometter, E. M. Spanton, J. I. A. Li, C. R. Dean, T. Taniguchi, K. Watanabe, M. Serbyn, and A. F. Young, Emergent Dirac gullies and gully-symmetry-breaking quantum Hall states in $ABA$ trilayer graphene, Phys. Rev. Lett. 121, 167601 (2018).
  * Note [1] This approach is aptly referred to as graduate-student descent.
  * Carleo  _et al._ [2019] G. Carleo, I. Cirac, K. Cranmer, L. Daudet, M. Schuld, N. Tishby, L. Vogt-Maranto, and L. Zdeborová, Machine learning and the physical sciences, Reviews of Modern Physics 91, 045002 (2019).
  * Carleo and Troyer [2017] G. Carleo and M. Troyer, Solving the quantum many-body problem with artificial neural networks, Science 355, 602 (2017).
  * Choo  _et al._ [2018] K. Choo, G. Carleo, N. Regnault, and T. Neupert, Symmetries and many-body excitations with neural-network quantum states, Physical review letters 121, 167204 (2018).
  * Cai and Liu [2018] Z. Cai and J. Liu, Approximating quantum many-body wave functions using artificial neural networks, Phys. Rev. B 97, 035116 (2018).
  * Borin and Abanin [2020] A. Borin and D. A. Abanin, Approximating power of machine-learning ansatz for quantum many-body states, Phys. Rev. B 101, 195141 (2020).
  * Christensen  _et al._ [2020] T. Christensen, C. Loh, S. Picek, D. Jakobovic, L. Jing, S. Fisher, V. Ceperic, J. D. Joannopoulos, and M. Soljacic, Predictive and generative machine learning models for photonic crystals, Nanophotonics 9, 4183 (2020).
  * Nikulin  _et al._ [2022] A. Nikulin, I. Zisman, M. Eich, A. Y. Petrov, and A. Itin, Machine learning models for photonic crystals band diagram prediction and gap optimisation, Photonics and Nanostructures - Fundamentals and Applications 52, 101076 (2022).
  * Sheverdin  _et al._ [2020] A. Sheverdin, F. Monticone, and C. Valagiannopoulos, Photonic inverse design with neural networks: The case of invisibility in the visible, Phys. Rev. Appl. 14, 024054 (2020).
  * Rem  _et al._ [2019] B. S. Rem, N. Käming, M. Tarnowski, L. Asteria, N. Fläschner, C. Becker, K. Sengstock, and C. Weitenberg, Identifying quantum phase transitions using artificial neural networks on experimental data, Nat. Phys. 15, 917 (2019).
  * Zhang  _et al._ [2019] Y. Zhang, A. Mesaros, K. Fujita, S. Edkins, M. Hamidian, K. Ch’ng, H. Eisaki, S. Uchida, J. S. Davis, E. Khatami, and E.-A. Kim, Machine learning in electronic-quantum-matter imaging experiments, Nature 570, 484 (2019).
  * Goodfellow  _et al._ [2016] I. Goodfellow, Y. Bengio, and A. Courville,  _Deep Learning_ (MIT Press, 2016) http://www.deeplearningbook.org.
  * Bishop [2006] C. M. Bishop,  _Pattern Recognition and Machine Learning_ (Springer-Verlag, 2006).
  * Hornik  _et al._ [1989] K. Hornik, M. Stinchcombe, and H. White, Multilayer feedforward networks are universal approximators, Neural Networks 2, 359 (1989).
  * Kidger and Lyons [2020] P. Kidger and T. Lyons, Universal Approximation with Deep Narrow Networks, in  _Proceedings of Thirty Third Conference on Learning Theory_ , Vol. 125, edited by J. Abernethy and S. Agarwal (PMLR, 2020) pp. 2306–2327.
  * Vapnik [1991] V. Vapnik, Principles of risk minimization for learning theory, Advances in neural information processing systems 4 (1991).
  * Dresselhaus and Dresselhaus [1981] M. Dresselhaus and G. Dresselhaus, Intercalation compounds of graphite, Advances in Physics 30, 139 (1981).
  * Serbyn and Abanin [2013] M. Serbyn and D. A. Abanin, New Dirac points and multiple Landau level crossings in biased trilayer graphene, Phys. Rev. B 87, 115422 (2013).
  * Glorot and Bengio [2010] X. Glorot and Y. Bengio, Understanding the difficulty of training deep feedforward neural networks, in  _Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics_, Proceedings of Machine Learning Research, Vol. 9, edited by Y. W. Teh and M. Titterington (PMLR, 2010) pp. 249–256.
  * Kingma and Ba [2015] D. P. Kingma and J. L. Ba, Adam: A method for stochastic optimization, in  _International Conference on Learning Representations (ICLR)_ (2015).
  * Clevert  _et al._ [2015] D.-A. Clevert, T. Unterthiner, and S. Hochreiter, Fast and accurate deep network learning by exponential linear units (ELUs) (2015), arXiv:1511.07289 .
  * Ba  _et al._ [2016] J. L. Ba, J. R. Kiros, and G. E. Hinton, Layer normalization (2016), arXiv:1607.06450 .
  * He  _et al._ [2016] K. He, X. Zhang, S. Ren, and J. Sun, Deep residual learning for image recognition, in  _2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ (2016) pp. 770–778.
  * Tancik  _et al._ [2020] M. Tancik, P. P. Srinivasan, B. Mildenhall, S. Fridovich-Keil, N. Raghavan, U. Singhal, R. Ramamoorthi, J. T. Barron, and R. Ng, Fourier features let networks learn high frequency functions in low dimensional domains, in  _NeurIPS_ (2020).
  * Abadi  _et al._ [2015] M. Abadi, A. Agarwal, P. Barham, E. Brevdo, Z. Chen, C. Citro, G. S. Corrado, A. Davis, J. Dean, M. Devin, S. Ghemawat, I. Goodfellow, A. Harp, G. Irving, M. Isard, Y. Jia, R. Jozefowicz, L. Kaiser, M. Kudlur, J. Levenberg, D. Mané, R. Monga, S. Moore, D. Murray, C. Olah, M. Schuster, J. Shlens, B. Steiner, I. Sutskever, K. Talwar, P. Tucker, V. Vanhoucke, V. Vasudevan, F. Viégas, O. Vinyals, P. Warden, M. Wattenberg, M. Wicke, Y. Yu, and X. Zheng, TensorFlow: Large-scale machine learning on heterogeneous systems (2015).
  * [37] https://git.ista.ac.at/qdyn/ML-band-structure-DOS
  * Note [2] We also experimented with the constrained trust-region method of Ref. [52], and the quasi-Newton method L-BFGS-B [53], both as implemented in scipy [54]. However, noise in the experimental data caused both these methods to perform poorly, becoming trapped in local optima.
  * Rumelhart  _et al._ [1986] D. E. Rumelhart, G. E. Hinton, and R. J. Williams, Learning representations by back-propagating errors, Nature 323, 533 (1986).
  * Hunt  _et al._ [2017] B. Hunt, J. Li, A. Zibrov, L. Wang, T. Taniguchi, K. Watanabe, J. Hone, C. Dean, M. Zaletel, R. Ashoori, and A. Young, Direct measurement of discrete valley and orbital quantum numbers in bilayer graphene, Nature comm. 8, 1 (2017).
  * Zibrov  _et al._ [2017] A. A. Zibrov, C. Kometter, H. Zhou, E. Spanton, T. Taniguchi, K. Watanabe, M. Zaletel, and A. Young, Tunable interacting composite fermion phases in a half-filled bilayer-graphene Landau level, Nature 549, 360 (2017).
  * Note [3] Note that the validation set is used for determining convergence and selecting the network architecture, while all results are reported on the testing set; this avoids overfitting to the testing data.
  * Haralick and Shapiro [1992] R. M. Haralick and L. G. Shapiro,  _Computer and Robot Vision, Volume II_ (Addison-Wesley, 1992) pp. 316–317.
  * Lewis [1995] J. P. Lewis, Fast template matching, Vision Interface , 120 (1995).
  * Studholme  _et al._ [1996] C. Studholme, D. L. G. Hill, and D. J. Hawkes, Automated 3-D registration of MR and CT images of the head, Medical Image Analysis 1, 163 (1996).
  * Klein*  _et al._ [2010] S. Klein*, M. Staring*, K. Murphy, M. A. Viergever, and J. P. Pluim, elastix: a toolbox for intensity-based medical image registration, IEEE Transactions on Medical Imaging 29, 196 (2010).
  * Brown [1992] L. G. Brown, A survey of image registration techniques, ACM Computing Surveys 24, 325 (1992).
  * Anderson  _et al._ [1984] C. H. Anderson, J. R. Bergen, P. J. Burt, and J. M. Ogden, Pyramid methods in image processing, RCA Engineer 29, 33 (1984).
  * Zhou  _et al._ [2021] H. Zhou, T. Xie, A. Ghazaryan, T. Holder, J. R. Ehrets, E. M. Spanton, T. Taniguchi, K. Watanabe, E. Berg, M. Serbyn, and A. F. Young, Half-and quarter-metals in rhombohedral trilayer graphene, Nature 598, 429 (2021).
  * Zondiner  _et al._ [2020] U. Zondiner, A. Rozen, D. Rodan-Legrain, Y. Cao, R. Queiroz, T. Taniguchi, K. Watanabe, Y. Oreg, F. von Oppen, A. Stern, E. Berg, P. Jarillo-Herrero, and S. Ilani, Cascade of phase transitions and Dirac revivals in magic-angle graphene, Nature 582, 203 (2020).
  * Wu  _et al._ [2021] S. Wu, Z. Zhang, K. Watanabe, T. Taniguchi, and E. Y. Andrei, Chern insulators, van Hove singularities and topological flat bands in magic-angle twisted bilayer graphene, Nature Mat. 20, 488 (2021).
  * Byrd  _et al._ [1999] R. H. Byrd, M. E. Hribar, and J. Nocedal, An interior point algorithm for large-scale nonlinear programming, SIAM Journal on Optimization 9, 877 (1999).
  * Zhu  _et al._ [1997] C. Zhu, R. H. Byrd, and J. Nocedal, Algorithm 778: L-BFGS-B: FORTRAN routines for large scale bound constrained optimization, ACM Transactions on Mathematical Software 23, 550 (1997).
  * Virtanen  _et al._ [2020] P. Virtanen, R. Gommers, T. E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau, E. Burovski, P. Peterson, W. Weckesser, J. Bright, S. J. van der Walt, M. Brett, J. Wilson, K. J. Millman, N. Mayorov, A. R. J. Nelson, E. Jones, R. Kern, E. Larson, C. J. Carey, İ. Polat, Y. Feng, E. W. Moore, J. VanderPlas, D. Laxalde, J. Perktold, R. Cimrman, I. Henriksen, E. A. Quintero, C. R. Harris, A. M. Archibald, A. H. Ribeiro, F. Pedregosa, P. van Mulbregt, and SciPy 1.0 Contributors, SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python, Nature Methods 17, 261 (2020).