# Photoemission orbital tomography based on tight-binding approach: method and
application to $\pi$-conjugated molecules

Misa Nozaki<EMAIL_ADDRESS>Graduate School of Science and Engineering,
Chiba University, Yayoi-cho 1-33, Inage, Chiba 263-8522, Japan  Takehisa
Konishi<EMAIL_ADDRESS>Graduate School of Science and
Engineering, Chiba University, Yayoi-cho 1-33, Inage, Chiba 263-8522, Japan

###### Abstract

Conventional photoemission orbital tomography based on Fourier iterative
method enables us to extract a projected two-dimensional (2D) molecular
orbital from a 2D photoelectron momentum map (PMM) of planar $\pi$-conjugated
molecules in a single-orientation system, while not in a multi-orientation
system. In this work, we demonstrate photoemission orbital tomography for
$\pi$-conjugated molecules with a tight-binding ansatz (linear combination of
atomic orbitals). We analyze 2D PMMs of single-orientation pentacene/Ag(110)
and multi-orientation 3,4,9,10-perylenetetracarboxylic dianhydride/Ag(110) and
reproduce their three-dimensional highest occupied molecular orbitals. We
demonstrate that the PhaseLift algorithm can be used to analyze PMM including
experimental or theoretical uncertainties. With the 2D PMM for pentacene, we
simultaneously optimized the structure and the molecular orbital. The present
approach enables us to extract the three-dimensional orbitals and structures
of existing materials.

##  I Introduction

Angle-resolved photoemission spectroscopy (ARPES) is an indispensable
experimental method to obtain information on the electronic states. ARPES has
been intensively used to address the electronic states of organic thin films
adsorbed on solid surfaces [1, 2, 3, 4, 5]. The recently developed time- and
angle-resolved photoemission spectroscopy [6, 7] made it possible to record
real-time non-equilibrium electronic processes in organic thin films [8].

Puschnig and others [2, 9] developed a method to reproduce molecular orbitals
from photoelectron momentum map (PMM) by recovering the phase of the Fourier-
transformed orbital. This method is called photoemission orbital tomography
(POT). PMM is proportional to the absolute square of the Fourier-transformed
molecular orbital within the framework of the plane-wave final state
approximation and does not contain information on the Fourier phases. The
currently proposed POT [9, 10, 11, 12, 13] uses the Fourier iterative method
[14, 15] and its variants [16, 13] as phase retrieval algorithms. These POTs
have been mainly used to reconstruct two-dimensional (2D) projections of
molecular orbitals from 2D PMMs of single-orientated planar $\pi$-conjugated
molecules adsorbed on metallic substrates [9, 10, 11, 12, 13]. Lüftner et al.
also accomplished 3D reconstruction of frontier orbitals of single-orientation
3,4,9,10-perylenetetracarboxylic dianhydride (PTCDA) by varying the photon
energy [17]. One should, however, note that preparing such a 3D PMM requires
prudent calibration of the intensities of light for different energies [17,
18].

For microscopic insights into various organic semiconductors, further
developments are desired. First, reconstructing a reliable 3D molecular
orbital from a single set of 2D PMM is desired. From a 2D PMM, we can
reproduce at most a 2D projection of molecular orbital with the conventional
POT based on the Fourier iterative method. Kliuev et al. discussed the
possibility to reconstruct a full 3D orbital from a 2D PMM by using an ansatz,
but did not show concrete examples [12]. Second, even in epitaxial molecular
films, it is quite common to find not one but several symmetry-related
molecular orientations. Therefore, it is important to extend POT to multi-
orientation systems. The conventional POT based on the Fourier iterative
method cannot treat incoherent superposition of photoelectron intensities from
molecules with various orientations [19].

In this work, we propose a POT based on a tight-binding picture and apply it
to reproducing molecular orbitals of $\pi$-conjugated molecules on metal
substrate. We formulate the POT as a fitting of experimental PMM to
theoretical photoelectron intensity from molecular orbitals. By describing the
molecular orbitals with localized atomic orbitals centered on fixed atomic
positions, the problem becomes simpler. To solve the problem, we use the least
squares fitting and the PhaseLift algorithm [20, 21, 22]. We reproduce the 3D
highest occupied molecular orbitals (HOMOs) of single-orientation pentacene
and multi-orientation PTCDA [Fig. 1(a) and 1(c)] from their 2D PMMs. In the
case of pentacene, we also remove the constraint of fixed atomic positions by
simultaneously optimizing the molecular orbital coefficients and the molecular
structure with the experimental 2D PMM.

##  II Theory

###  II.1 Photoelectron intensity of molecules

In the plane wave final state approximation, we can describe the photoelectron
intensity from an adsorbed molecule on metal substrate as

$\displaystyle W_{{\bf k}}(\ket{\psi_{i}})$ $\displaystyle\propto$
$\displaystyle|\langle{\bf k}|{\bm{\varepsilon}}\cdot{\bf
p}|\psi_{i}\rangle|^{2}\delta\left(h\nu-E_{i}-E_{k}-\Phi\right)$ (1)
$\displaystyle\propto$ $\displaystyle|{\bm{\varepsilon}}\cdot{\bf
k}|^{2}|\braket{{\bf k}}{\psi_{i}}|^{2}\delta\left(h\nu-
E_{i}-E_{k}-\Phi\right),$

using Fermi’s golden rule. Here $W_{{\bf k}}(\ket{\psi_{i}})$ denotes the
transition rate from an initial state $\ket{\psi_{i}}$ with binding energy
$E_{i}$ relative to the Fermi energy of the substrate to a plane wave final
state $\ket{{\bf k}}$ with kinetic energy $E_{k}$. $h\nu$ and
${\bm{\varepsilon}}$ are the photon energy and the polarization vector of the
incident light, ${\bf p}$ is the momentum operator, $\Phi$ is the work
function. In this work, we focus on nondegenerate molecular orbitals and
consider real-valued $\psi_{i}({\bf r})$. When we discuss a 2D PMM for a fixed
$E_{k}=h\nu-E_{i}-\Phi$, the photoelectron intensity is proportional to

$\displaystyle I_{{\bf k}}(\ket{\psi_{i}})=|{\bm{\varepsilon}}\cdot{\bf
k}|^{2}|\braket{{\bf k}}{\psi_{i}}|^{2}.$ (2)

Let us apply Eq. (2) to an organic molecular layer. We consider systems
consisting of weakly interacting molecules in which the electronic states
retain the character of molecular orbitals of an isolated molecule. In such
systems, the photoelectrons emitted from individual molecules do not interfere
with each other, and the net photoelectron intensity is proportional to the
incoherent sum of $I_{{\bf k}}$’s in Eq. (2) for different molecules. Writing
the molecular orbital of the $r$-th molecule as $\ket{\psi_{i}^{r}}$ and
assuming the binding energy $E_{i}^{r}$ to be the same for all $r$
($r=1,2,\cdots,p$), the photoelectron intensity for the system is proportional
to

$\displaystyle I_{{\bf
k}}(\\{\ket{\psi_{i}^{r}}\\})=|{\bm{\varepsilon}}\cdot{\bf
k}|^{2}\sum_{r=1}^{p}|\braket{{\bf k}}{\psi_{i}^{r}}|^{2}.$ (3)

Here $\\{\ket{\psi_{i}^{r}}\\}$ is the set of orbitals for all the molecules.

We can describe Eq. (3) with a single molecular orbital when the orbitals of
different molecules are related to each other by translation and rotation.
Setting $\ket{\psi_{i}^{1}}$ as a reference, $\ket{\psi_{i}^{r}}$ is expressed
using translation $T_{r}$ and rotation $R_{r}$ as follows:

$\displaystyle\ket{\psi_{i}^{r}}=T_{r}R_{r}e^{i\theta_{r}}\ket{\psi_{i}^{1}}\quad(T_{1}=R_{1}=I).$
(4)

Here $I$ is the identity operator, and $e^{i\theta_{r}}$ is a phase factor.
Since translations of molecules do not affect the photoelectron intensity, we
drop $T_{r}$ below. Substituting Eq. (4) into Eq.(3), we obtain

$\displaystyle I_{{\bf k}}(\ket{\psi_{i}})=|{\bm{\varepsilon}}\cdot{\bf
k}|^{2}\sum_{R}w_{R}|\bra{{\bf k}}R\ket{\psi_{i}}|^{2}.$ (5)

Here $R$ indicates different operations, and $w_{R}$ the number of molecules
(or the number ratio of molecules) specified by $R$. We omit superscript 1 for
the reference system.

Finally, we express the photoelectron intensity (5) in terms of the projection
operator

$\displaystyle P=\ket{\psi_{i}}\bra{\psi_{i}},$ (6)

which is convenient for our POT. The intensity is

$\displaystyle I_{{\bf k}}(P)=|{\bm{\varepsilon}}\cdot{\bf k}|^{2}{\rm
Tr}\left[F_{{\bf k}}P\right],$ (7)

with

$\displaystyle F_{{\bf k}}=\sum_{R}w_{R}R^{-1}\ket{{\bf k}}\bra{{\bf k}}{R}.$
(8)

The projection operator ${P}$ fulfills

$\displaystyle\quad\langle\psi|{P}|\psi\rangle\geq
0,\quad\text{Tr}[P]=1,\quad\text{rank}({P})=1,$ (9)

where $|\psi\rangle$ is an arbitrary single electron state. We denote the
first property (positive semidefiniteness) as ${P}\succeq 0$.

To numerically handle Eq. (7), we approximate $\ket{\psi_{i}}$ as

$\displaystyle\ket{\psi_{i}}=\sum_{n=1}^{N}\ket{\chi_{n}}c_{n}$ (10)

with a finite real-valued orthonormal basis $\\{{\chi_{n}}({\bf
r})|n=1,...,N\\}$. Regarding the expansion coefficient $c_{n}$ as the $n$-th
component of a normalized real column vector ${\bf c}$,

$\displaystyle{{\bf
c}}=\left(c_{1},c_{2},\cdots,c_{n},\cdots,c_{N}\right)^{T},$ (11)

$P$ reduces to a matrix

$\displaystyle C={{\bf c}}{{\bf c}}^{T}.$ (12)

Thereby, we obtain the matrix representation of Eq. (7) in terms of $C$ as

$\displaystyle I_{{\bf k}}(C)$ $\displaystyle=$
$\displaystyle|{\bm{\varepsilon}}\cdot{\bf
k}|^{2}\sum_{n,n^{\prime}=1}^{N}\bra{\chi_{n}}F_{{\bf
k}}\ket{\chi_{n^{\prime}}}c_{n^{\prime}}c_{n}$ (13) $\displaystyle=$
$\displaystyle|{\bm{\varepsilon}}\cdot{\bf k}|^{2}\textrm{Tr}\left[F_{{\bf
k}}C\right],$

where $F_{{\bf k}}$ in the second line is a matrix whose
$n,n^{\prime}$-component is $\bra{\chi_{n}}F_{{\bf
k}}\ket{\chi_{n^{\prime}}}$. $C$ inherits the properties of $P$ in Eq. (9):

$\displaystyle C\succeq 0,\quad\text{Tr}[C]=1,\quad\text{rank}(C)=1.$ (14)

Eq. (13) with conditions Eq. (14) turns out to be a suitable expression for
the PhaseLift algorithm as shown below (Sec. II.3).

###  II.2 Photoemission orbital tomography based on least squares fitting

We aim at determining ${\bf c}$ for $\ket{\psi_{i}}$ by fitting experimental
photoelectron intensities to Eq. (13). Suppose that we have experimental
photoelectron intensities $z_{m}$ from an unknown molecular orbital
$\ket{\psi_{i}}$ at ${{\bf k}}={{\bf k}_{m}}$ $(m=1,2,\cdots,M)$. Our problem
is a non-linear least-squares problem,

$\displaystyle\underset{{\bf c},\lambda}{\text{Minimize}}\quad$
$\displaystyle\left\|{\bf z}-\lambda{\bf I}({\bf c}{\bf c}^{T})\right\|_{2}$
(15)

where ${\bf z}=(z_{1},z_{2},\cdots,z_{M})^{T}$ and ${\bf I}(A)=(I_{{\bf
k}_{1}}(A),I_{{\bf k}_{2}}(A)\cdots,I_{{\bf k}_{M}}(A))^{T}$. $\|\cdot\|_{2}$
represents the Euclidean distance ($L_{2}$ norm), and $\lambda$ is the
constant of proportionality.

To solve this problem, we use the trust region reflective (TRF) algorithm [23]
and the basin hopping algorithm [24]. In the least squares fitting, the
results can be at local minima depending on the initial values of $\bf c$.

###  II.3 Orbital tomography based on PhaseLift

For POT, we can also use a PhaseLift based algorithm. PhaseLift algorithm is
classified as a semidefinite program (SDP) and is known to be robust against
additive noise [20].

In the PhaseLift-based approach, we find a low-rank real symmetric matrix
$\overline{C}$ instead of $\lambda$ and $\bf c$ by minimizing
$\text{Tr}[{\overline{C}}]$ under the experimental constraints. We solve the
following problem with a finite uncertainty $\eta$,

$\displaystyle\underset{\overline{C}}{\text{Minimize}}\quad$
$\displaystyle{\text{Tr}[\overline{C}]}$ subject to $\displaystyle\left\|{{\bf
z}}-{{\bf I}}(\overline{C})\right\|_{2}\leq\eta,$
$\displaystyle\overline{C}\succeq 0.$ (16)

In reality, there can be contributions to $\eta$ other than the statistical
noise, such as systematic uncertainties in the data due to, for example,
signals from the substrate, or radiation damage of the sample. There are also
uncertainties which arise from the theoretical model such as use of a
insuficient basis function set, errors in the assumed molecular structure, or
the plane wave approximation of the photoelectron final states.

The rank of resulting $\overline{C}$ can be larger than 1. In that case,
$\overline{C}$ represents a mixed state containing more than one component.
Following Ref. [20], we extract the rank 1 matrix $C$ by eigenvalue
decomposition

$\displaystyle\overline{C}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\nu_{i}{\bf u}_{i}{\bf
u}_{i}^{T},\quad\nu_{1}\geq\nu_{2}\geq\cdots\geq\nu_{N}\geq 0.$ (17)

Here $\nu_{i}$ and ${\bf u}_{i}$ are the eigenvalues and normalized
eigenvectors of the matrix $\overline{C}$. We regard ${\bf u}_{1}$ as $\bf c$,

$\displaystyle C={\bf u}_{1}{\bf u}_{1}^{T}.$ (18)

In order to correct for the loss of the intensity caused by the disregard of
the minor components, we rescale the intensity by minimizing the distance
between $\bf z$ and rescaled ${\bf I}(C)$, $\left\|{\bf z}-\lambda{\bf
I}(C)\right\|_{2}$. The distance is minimum for $\lambda=\frac{{\bf
z}\cdot{\bf I}(C)}{\left\|{\bf I}(C)\right\|_{2}^{2}}$ where it has the value

$\displaystyle d(C)={\left\|{\bf z}\right\|_{2}}\left(1-\frac{\left({\bf
z}\cdot{\bf I}(C)\right)^{2}}{{\left\|{\bf z}\right\|_{2}^{2}}\left\|{\bf
I}(C)\right\|_{2}^{2}}\right)^{1/2}.$ (19)

Replacement of $\nu_{1}$ by the optimized $\lambda$ is similar to the
“debiasing” of the solution in Ref. [20]. This is also useful when
$\overline{C}$ is rank 1 because intensities tend to be scaled down through
trace minimization, especially at large $\eta$. Since the appropriate value of
$\eta$ for each particular case is not known a priori, as a practical means
for choosing $\eta$, we find $C$ with the smallest $d(C)$ by repeating the
same procedure with different $\eta$’s.

###  II.4 Atomic orbital basis set

Since the molecular orbitals of our interest are localized around the
molecule, they can be approximated with a small number of atomic orbitals
[25]. Atomic orbitals $\\{\ket{\phi_{n}}|n=1,2,\cdots,N\\}$ are usually
nonorthogonal, yet we can easily orthonormalize them. To this end, we can use,
for example, Löwdin’s linear transformation [26, 25]. Diagonalizing the
overlap matrix $S_{nn^{\prime}}=\langle\phi_{n}|\phi_{n^{\prime}}\rangle$,

$\displaystyle U^{T}SU$ $\displaystyle=s,$ (20)

we define matrix V as

$\displaystyle V$ $\displaystyle=Us^{-\frac{1}{2}}.$ (21)

$U$ is an orthogonal matrix composed of the eigenvectors of $S$, $s$ the
diagonal matrix whose diagonal elements are the eigenvalues of $S$, and
$s^{-\frac{1}{2}}$ is the diagonal matrix of the inverse square root of the
eigenvalues of $S$. With $V$, we can orthonormalize the atomic orbitals as

$\displaystyle\ket{\phi_{n}^{\prime}}=\sum_{q}\ket{\phi_{q}}V_{qn}.$ (22)

We use these orthonormal orbitals $\\{\ket{\phi_{n}^{\prime}}\\}$ as
$\\{\ket{\chi_{n}}\\}$ in Eq. (10).

##  III Application

In this work, we apply the above method to the highest occupied molecular
orbital (HOMO) of monolayer pentacene on Ag(110) (single-orientation system)
and the HOMO of herringbone monolayer PTCDA on Ag(110) surfaces (multi-
orientation system). Figures 1(a) and 1(c) illustrate surface models of
monolayer pentacene/Ag(110) and monolayer PTCDA/Ag(110), respectively. All the
pentacene molecules lie flat on the silver surface and are aligned with their
long molecular axis along the [001] direction of the silver substrate [27, 4].
All the PTCDA molecules are lying flat on the silver surface, half of them
aligned with their long molecular axis along the [001] and the other half
along the [1$\bar{1}$0] direction of the silver substrate [28, 29]. To our
knowledge, their precise adsorption positions have not been determined
experimentally. In both systems, the HOMO does not show significant
hybridization with the substrate [30, 31, 29, 4].

###  III.1 Experimental data

####  III.1.1 pentacene/Ag(110)

As the experimental PMM of the HOMO of pentacene/Ag(110), we take the
experimental data set by Metzger et al. reported in Ref. [32]. In their
experiment, photoelectrons excited from the initial state at the binding
energy $E_{i}$ = 1.2 eV with the circularly polarized light of $h\nu$ = 28 eV
incident from the polar angle $\theta$ = 65∘ and azimuthal angle $\phi$ = -85∘
are collected [Fig. 1(b)]. As the value of the work function is not reported
in Ref. [32], we use the value $\Phi$ = 4.0 eV and $E_{k}$ = 22.8 eV. This
arbitrary choice of $E_{k}$ does not have a major effect on the orbital
reconstruction except for the quantitative molecular structure determination
described in Sec.IV.3. The data in Ref. [32] have been background-substracted
and symmetrized with respect to $k_{x}$ and $k_{y}$ axes.

We prepared the input experimental PMM data for our analysis shown in Fig.
3(a) from the data in Ref. [32] in two steps. In the first step, we extracted
numerical values from the 2D color map of PMM [Fig. 4.5(b) in Ref. [32]]. In
the second step, we performed an average pooling over the range of $\Delta
k_{x}\Delta k_{y}$ $\approx$ 0.072 Å-2 on the digitized experimental PMM.

In the estimation of the molecular orbital, it is necessary to take into
account the experimental condition and the effect of the above-mentioned data
processing. The symmetrization of the data is equivalent to the superposition
of the four PMMs with incident light from
$(\theta,\phi)=(65^{\circ},-85^{\circ}),(65^{\circ},85^{\circ}),(65^{\circ},95^{\circ})$,
and $(65^{\circ}-95^{\circ})$. We take the sum of
$|{\bm{\varepsilon}}\cdot{\bf k}|^{2}$ for the four incident directions in Eq.
(13).

####  III.1.2 PTCDA/Ag(110)

As the experimental PMM for the HOMO of PTCDA/Ag(110), we refer to the
experimental data by Willenbockel et al. reported in Ref. [28]. In their
experiment, the photoelectrons are excited from the initial state at $E_{i}$
=1.74 eV by the $p$-polarized light of $h\nu$ = 35 eV with incident angle
$\theta$ = 40∘. The sample was rotated around the $z$-axis, and the
photoelectrons emitted in the plane spanned by the $z$-axis and the direction
of the incident light were detected [Fig. 1(d)]. As the value of $\Phi$ is not
reported in Ref. [28], we use the value $\Phi$ = 4.0 eV and $E_{k}$ = 29.3 eV.

We digitized the PMM in the same manner as pentacene/Ag(110) described in Sec.
III.1.1. In our orbital estimations, we assume that the orientations of the
two molecules in the surface unit cell are related to each other by 90∘
rotation. We assume the PMM is symmetric with respect to the line $k_{x}$ =
$k_{y}$ and symmetrized the digitized data accordingly [Fig. 5(a)].

###  III.2 Computational details

####  III.2.1 Density functional theory calculation

We performed structural optimization for isolated pentacene and PTCDA at the
level of B3LYP [33] with the 6-31++G(d,p) basis set using Gaussian16 [34]. The
obtained structures were used in the estimations of molecular orbitals (Sec.
IV.1.1 and Sec. IV.1.2), and the Kohn-Sham (KS) orbitals were compared with
the estimated orbitals [Fig. 3(e) and Fig. 5(d)].

For comparison with the optimized structure using PMM, we also performed the
structural optimization for pentacene monolayers adsorbed on Ag(110) surface
using the Vienna Ab initio Simulation Package (VASP) [35]. We used the
projector-augmented wave (PAW) [36] approach with an energy cutoff of 400 eV.
The PBE functional [37] was chosen for exchange-correlation effects, and van
der Waals corrections are treated with the DFT-D3 method of Grimme with zero-
damping function [38].

We performed the calculations for two possible adsorption positions [Fig. 2].
The surface structure was modeled using a 10-layer Ag(110) slab with a
pentacene molecule separated by a vacuum layer of 15 Å. The 2D lattice vectors
are described using the superstructure matrix $\begin{pmatrix}3&-1\\\
1&4\end{pmatrix}$ [3] with the substrate 2D lattice vectors (a,0) and
(0,$\sqrt{2}$a) (a = 4.07 Å).

In the electronic structure calculations, we used a first-order Methfessel-
Paxton smearing [39] of 0.2 eV, 3$\times$3$\times$1 Monkhorst-Pack k-mesh grid
[40], and the energy convergence tolerance of 10-7 eV/unit cell. We considered
the electronic dipole correction along the c-axis to take into account the
possible effect of charge interactions between the repeated slabs. In the
structural optimizations, all the ionic positions were optimized until all
residual forces were less than 10-2 eV/Å.

(a) | (b)
---|---
|
(c) | (d)
|
Figure 1: (a),(c) Structural models of (a) pentacene/Ag(110) and (c) PTCDA/Ag(110). The dark gray, black, white, and light gray spheres are, respectively, carbon, oxygen, hydrogen, and silver atoms. The parallelograms in (a) and (c) represent surface unit cells. (b),(d) The relation between the molecular coordinates and the incident direction of the light for (b) pentacene/Ag(110) and (d) PTCDA/Ag(110) in the measurement in Ref. [32] and Ref. [28]. For PTCDA/Ag(110) system, the incident direction of the light and photoelectron wave vector ${\bf k}$ reside in the same plane perpendicular to the molecular plane. (a) | (b)
---|---
|
Figure 2: Model surface structures of pentacene/Ag(110). In (a) and (b), the
centroid of pentacene is at the hollow site and short-bridge site,
respectively.

####  III.2.2 POT

As the input molecular structures, we used the optimized structures of
isolated pentacene and PTCDA within the density functional theory (DFT). As
the atomic orbitals $\phi_{q}$ in Eq. (22), we used the STO-3G basis set and
considered only 2$p_{z}$ orbitals of carbon and oxygen atoms in the molecules,
omitting the other basis functions in the set. The ground for using the
minimal basis set is given in Appendix. In the orbital estimation of PTCDA
with herringbone structure, $R$ in Eq. (8) includes the identity operator and
the 90∘ rotation operator around the z-axis. $w_{R}$ for both operators are
fixed to be 1. In the least squares fitting, we checked that the TRF and the
basin hopping algorithms give the same result by adequately choosing the
initial values. In the PhaseLift based approach, we varied $\eta^{\prime}$
($\eta=\eta^{\prime}\cdot\|{\bf z}\|_{2}$) from 0 to 0.95 by 0.05 step and
solved problem (16) using Splitting Conic Solver (version: 3.1.0) [41] via
cvxpy (version: 1.1.17) [42, 43].

In the orbital reconstruction, we imposed the constraints that the molecular
orbital coefficients are symmetric or anti-symmetric with respect to the
reflection through the $yz$ and $xz$ plane. These constraints, which are
required by the structural symmetry of both pentacene and PTCDA, were
necessary for the orbital reconstruction from the present data set.

##  IV Results and Discussion

###  IV.1 Orbital estimation

####  IV.1.1 Pentacene/Ag(110)

We estimated the HOMO of pentacene by using the least squares approach (Sec.
II.2). Figures 3(b) and 3(f) show the fitted PMM and the estimated orbital. We
could well reproduce the experimental PMM with the estimated orbital. The
normalized distance between the experimental PMM [Fig. 3(a)] and the fitted
PMM $\|{\bf z}-\lambda{\bf I}({\bf c}{\bf c}^{T})\|_{2}/\|{\bf z}\|_{2}$ is
0.171. The spatial pattern of the estimated orbital and the KS-HOMO [Fig.
3(e)] are consistent with each other.

We also estimated the HOMO of pentacene by using the PhaseLift based approach
(Sec. II.3). $d(C)/\|{\bf z}\|_{2}$ takes the minimum value 0.175 at
$\eta^{\prime}$ = 0.25 (Fig. 4). Figures 3(c) and 3(g) show the fitted PMM
$\lambda{\bf I}(C)$ and the estimated orbital at $\eta^{\prime}$ = 0.25. We
could well reproduce the experimental PMM with the estimated orbital [Fig.
3(g)]. The estimated orbital resembles KS-HOMO [Fig. 3 (e)] as well as the
estimated orbital of the least squares approach.

(a) | (b) | (c) | (d)
---|---|---|---
|  |  |
(e) | (f) | (g) | (h)
|  |  |
Figure 3:  (a) Experimental PMM of the HOMO of pentacene/Ag(110) (original
data were taken from Fig. 4.5(b) in Ref. [32]). (b) The fitted PMM with the
least squares approach. (c) The fitted PMM $\lambda{\bf I}(C)$ with the
PhaseLift approach at $\eta^{\prime}$ = 0.25. In (b) and (c), the DFT-
optimized molecular structure was used. (d) The fitted PMM with the least
squares approach for the PMM-optimized structure. The plotted data are
normalized as ${\bf z}/\|{\bf z}\|_{2}$ for (a) and $\lambda{\bf I}(C)/\|{\bf
z}\|_{2}$ for (b), (c) and (d). (e) The KS-HOMO of an isolated pentacene. (f)
The estimated orbital with the least squares approach. (g) The estimated
orbital with the PhaseLift based approach at $\eta^{\prime}$ = 0.25. (h) The
estimated orbital with the least squares approach for the molecular structure
optimized using PMM. In (e)-(h), the isosurface value is 0.02 a.u. The
orbitals are visualized using the VESTA software [44].  Figure 4: $d(C)$ and
eigenvalues of $\overline{C}$ for pentacene/Ag(110). Crosses show the distance
between the experimental PMM and fitted PMM ($d(C)$, left axis). Open symbols
show the largest three eigenvalues of $\overline{C}$ (right axis). In the
shaded region, we could not obtain $\overline{C}$.

####  IV.1.2 PTCDA/Ag(110)

In the previous application on pentacene, we have shown that with our method
we can obtain the three-dimensional molecular orbital in a singe-orientational
system, with an accuracy comparable to the Fourier-based POT. In this section,
we apply our method to a multi-orientational system, namely PTCDA on Ag(110),
where the Fourier-based POT cannot be applied.

Figures 5(b) and 5(e) show the fitted PMM and the estimated orbital obtained
from the least squares approach. With the estimated orbital, the corner-frame-
like features of the experimental PMM [Fig. 5(a)] is well reproduced, whereas
the features around the $k_{x}$ = $k_{y}$ = 0 Å-1 are not reproduced. The
normalized distance between the experimental and the fitted PMMs $\|{\bf
z}-\lambda{\bf I}({\bf c}{\bf c}^{T})\|_{2}/\|{\bf z}\|_{2}$ is 0.379. Some of
the node positions of the estimated orbital [Fig. 5(e)] and the KS-HOMO [Fig.
5(d)] are different.

In the PhaseLift based approach, $d(C)/\|{\bf z}\|_{2}$ takes the minimum
value 0.486 at $\eta^{\prime}$ = 0.45 [Fig. 6(a)]. Figures 5(c) and 5(f) show
the fitted PMM $\lambda{\bf I}(C)$ and the estimated orbital at
$\eta^{\prime}$ = 0.45. With the estimated orbital, the corner-frame-like
features of the experimental PMM is well reproduced. The estimated orbital is
qualitatively different from the least-squares result, but resembles the KS-
HOMO [Fig. 5(d)] and reconstructed HOMO in Ref. [17]. The origin of the
features that were not reproduced by the principal component of $\overline{C}$
is unclear. The second and the third components at $\eta^{\prime}$ = 0.45
[shown in Fig. 6(b) and 6(c)] do not coincide with any of the Kohn-Sham
orbitals.

(a) | (b) | (c)
---|---|---
|  |
(d) | (e) | (f)
|  |
Figure 5:  (a) Experimental PMM of the HOMO of PTCDA/Ag(110) (original data ware taken from Fig. 2(c) in Ref. [28]). (b) The fitted PMM with the least squares approach. (c) The fitted PMM $\lambda{\bf I}(C)$ with the PhaseLift approach at $\eta^{\prime}$ = 0.25. In (b) and (c), the DFT-optimized structure was used. The plotted data are normalized as ${\bf z}/\|{\bf z}\|_{2}$ for (a) and $\lambda{\bf I}(C)/\|{\bf z}\|_{2}$ for (b) and (c). (d) The KS-HOMO of isolated PTCDA. (e) The estimated orbital with the least squares approach. (f) The estimated orbital with PhaseLift approach at $\eta^{\prime}$ = 0.45. In (d)-(f), the isosurface value is 0.02 a.u.  (a) |
---|---
(b) | (c)
|
Figure 6:  (a) $d(C)$ and eigenvalues of $\overline{C}$ for PTCDA/Ag(110).
Crosses show the distance between the experimental PMM and fitted PMM ($d(C)$,
left axis), while open symbols show the largest three eigenvalues of
$\overline{C}$ (right axis). In the shaded region, we could not obtain
$\overline{C}$. (b) and (c) represent components with the second and third
largest eigenvalues at $\eta^{\prime}$ = 0.45, respectively.

###  IV.2 Simultaneous optimization of molecular structure

One drawback of the present tight-binding based method is that the precise
molecular structure is necessary as an input, which is not always at hand. In
this regard, it is possible to remove the constraints of fixed atomic
positions. We performed the simultaneous molecular orbital/structural
optimization of pentacene using the experimental PMM of pentacene/Ag(110)
[Fig. 3(a)]. We estimated the molecular orbital coefficients using the TRF
algorithm and updated the structure using the steepest descent method. We
stopped the optimization when the difference of $\left\|{\bf z}-\lambda{\bf
I}({\bf c}{\bf c}^{T})\right\|_{2}/\|{\bf z}\|_{2}$ is less than $10^{-5}$. We
considered possible deformations of the molecule that do not change the D2h
symmetry.

After the PMM-based structural optimization, the normalized distance between
the fit and experimental PMM decreased from 0.171 to 0.166. Figure 3(d) shows
the fitted PMM with the optimized structure. The spatial pattern of the
obtained orbital [Fig. 3(h)] does not change much from Figs. 3(e), (f) and
(g). Thus, we can optimize atomic positions as well as molecular orbital
coefficients utilizing the information contained in the PMM.

###  IV.3 POT as molecular structure analysis

Given the above results, it is tempting to assess the possibility to extract
quantitative molecular structural information from PMM data using the present
tight-binding based POT method. This can be another advantage of using
localized atomic orbital basis functions, i.e., each of the atomic orbitals
used in the fitting has well-defined center at the atomic position in the
molecule and therefore we can expect that the quality of the fit is sensitive
to the atomic positions in the structural model.

Since PMM is essentially a Fourier transform of molecular orbital, accurate
absolute value of the photoelectron wavenumber is necessary to set the
absolute length scale for structural analysis. For this purpose, it may be
possible to improve the accuracy of the result by introducing a correction to
the kinetic energy measured from the vacuum level, that is, if the
photoemission from a $\pi$-conjugated molecule is approximately viewed as an
escape from a constant potential well, it would be necessary to make a
correction to $E_{k}$ by shifting the origin of the kinetic energy
accordingly, corresponding to the inner potential usually introduced as a
parameter in the case of solid surfaces.

In order to check the dependence of the result of the structural optimization
on the shift of $E_{k}\rightarrow E_{k}+\Delta E_{k}$, we performed structural
optimization with different values of $\Delta E_{k}$ (= -1 eV, 0 eV, 1 eV, 2
eV). The calculation with $\Delta E_{k}$ = 0 eV is equivalent to the
calculation in Sec. IV.2. The resulting distances C1-C2 and C3-C4 (defined in
Fig.2) from PMM-optimized structures are shown in Table.1. As expected from
the property of Fourier transformation, the resulting distances decrease as
$\Delta E_{k}$ increases. On the other hand, we obtained almost the same
$\left\|{\bf z}-\lambda{\bf I}({\bf c}{\bf c}^{T})\right\|_{2}/\|{\bf
z}\|_{2}$ with these values of $\Delta E_{k}$. We need additional data, e.g.,
$h\nu$ dependence of the intensity in order to determine $\Delta E_{k}$.

The results themselves are promising. From the $\Delta E_{k}$ dependence, the
deviation of the distances corresponding to an uncertainty of $\sim$1 eV in
the kinetic energy is $\sim$0.1% in the direction of the long molecular axis
and $\sim$3% in the direction of the short molecular axis. Also, we can see
that the results of the DFT calculation on pentacene/Ag(110) and the
experimental value for crystal fall in the range of the values resulting from
the PMM optimization, indicating that the PMM optimization gives values in the
plausible range and the PMM data may well contain relevant structural
information.

Table 1: PMM-optimized (PMM-opt) and DFT-optimized (DFT-opt) structures of pentacene/Ag(110), and experimental (Exp) structure of crystalline pentacene. For the DFT calculations, we used the configurations in Fig.2(a) and 2(b). The experimentally determined crystal structure in fiber-structured thin films is taken from Ref. [45]. method |  $\Delta E_{k}$ (eV) | C1-C2 (Å) | C3-C4 (Å) | $\frac{\|\lambda{\bf I}({\bf c}{\bf c}^{T})-{\bf z}\|_{2}}{\|{\bf z}\|_{2}}$
---|---|---|---|---
PMM-opt | -1.0 | 12.33 | 3.03 | 0.166
PMM-opt | 0.0 | 12.29 | 2.93 | 0.165
PMM-opt | 1.0 | 12.27 | 2.84 | 0.166
PMM-opt | 2.0 | 12.24 | 2.74 | 0.166
DFT-opt(a) |  | 12.33 | 2.81 |
DFT-opt(b) |  | 12.30 | 2.82 |
Exp. (crystal) |  | 12.19 | 2.81 |

##  V Conclusion

In summary, we have presented a photoemission orbital tomography (POT) method
based on the LCAO ansatz for the molecular orbitals and showed that it
overcomes some of the limitations of the methods based on the Fourier
iterative reconstruction. We successfully reconstructed the 3D HOMOs of
pentacene and PTCDA on Ag(110) from 2D PMMs at a single photon energy. In the
case of PTCDA/Ag(110), we demonstrated that multi-orientational molecular
films can be handled, where the PMMs of more than one molecular orientations
are superimposed. The comparison between the PhaseLift based algorithm and the
least squares algorithms showed that the former gives a more plausible orbital
in the case of PTCDA/Ag(110). Furthermore, we showed that the method can be
used to refine the molecular structure of the pentacene with the experimental
PMM.

In principle, information on intermolecular and molecule-substrate orbital
hybridization may also be extracted from the PMMs by including orbitals from
the surroundings in the basis set. The present POT method should become a
valuable tool for the microscopic understanding of the electronic interactions
governing the electronic properties of organic semiconductors.

###### Acknowledgements.

We would like to thank Manabu Hagiwara, Kazushi Mimura, and Kaori Niki for the
fruitful discussions. We in particular thank Kazushi Mimura for providing us
with a sample code to use the PhaseLift algorithm via cvxpy. M.N. thanks Peter
Krüger and Naoya Iwahara for careful reading of the manuscript. This work was
supported by Japan Science and Technology Agency, the Establishment of
University Fellowships Towards the Creation of Science Technology Innovation,
Grant Number JPMJFS2107.

*

##  Appendix A Use of larger basis sets in the PMM fitting

In general, it is not possible to improve the orbital reconstruction from 2D
PMM by use of split-valence basis sets such as multiple-zeta basis sets.
Suppose that we represent a $\pi$-conjugated molecular orbital as

$\psi_{i}({\bf r})=\sum_{n}Y_{10}(\widehat{{\bf r}-{\bf
R}_{n}})\sum_{\nu=1}^{N_{p}}d_{n\nu}R_{1\nu}(|{\bf r}-{\bf R}_{n}|),$ (23)

using $N_{p}$ radial functions for the $2p_{z}$ orbital. Its Fourier transform
is

${\tilde{\psi}}_{i}({\bf k})=-4\pi iY_{10}(\hat{\bf k})\sum_{n}e^{-i{\bf
k}\cdot{\bf R}_{n}}\sum_{\nu=1}^{N_{p}}d_{n\nu}\tilde{R}_{1\nu}(k),$ (24)

where

$\tilde{R}_{1\nu}(k)=\int_{0}^{\infty}drr^{2}j_{1}(kr)R_{1\nu}(r).$ (25)

It can be seen from Eq. (24) that it is impossible to determine the
coefficients $d_{n\nu}$’s at each atomic site, i.e., for the same $n$ and
different $\nu$ independently from a data set with fixed $|{\bf k}|=k$.
Optimizing $d_{n\nu}$ will be possible with 2D PMMs obtained at several photon
energies.

## References

  * Nakayama  _et al._ [2020] Y. Nakayama, S. Kera, and N. Ueno, Photoelectron spectroscopy on single crystals of organic semiconductors: experimental electronic band structure for optoelectronic properties, J. Mater. Chem. C 8, 9090 (2020).
  * Puschnig  _et al._ [2009] P. Puschnig, S. Berkebile, A. J. Fleming, G. Koller, K. Emtsev, T. Seyller, J. D. Riley, C. Ambrosch-Draxl, F. P. Netzer, and M. G. Ramsey, Reconstruction of Molecular Orbital Densities from Photoemission Data, Science 326, 702 (2009).
  * Puschnig  _et al._ [2011] P. Puschnig, E.-M. Reinisch, T. Ules, G. Koller, S. Soubatch, M. Ostler, L. Romaner, F. S. Tautz, C. Ambrosch-Draxl, and M. G. Ramsey, Orbital tomography: Deconvoluting photoemission spectra of organic molecules, Phys. Rev. B 84, 235427 (2011).
  * Ules  _et al._ [2014] T. Ules, D. Lüftner, E. M. Reinisch, G. Koller, P. Puschnig, and M. G. Ramsey, Orbital tomography of hybridized and dispersing molecular overlayers, Phys. Rev. B 90, 155430 (2014).
  * Iwasawa  _et al._ [2022] M. Iwasawa, S. Kobayashi, M. Sasaki, Y. Hasegawa, H. Ishii, F. Matsui, S. Kera, and Y. Yamada, Photoemission Tomography of a One-Dimensional Row Structure of a Flat-Lying Picene Multilayer on Ag(110), J. Phys. Chem. Lett. 13, 1512 (2022).
  * Kutnyakhov  _et al._ [2020] D. Kutnyakhov, R. P. Xian, M. Dendzik, M. Heber, F. Pressacco, S. Y. Agustsson, L. Wenthaus, H. Meyer, S. Gieschen, G. Mercurio, A. Benz, K. Bühlman, S. Ds̈ter, R. Gort, D. Curcio, K. Volckaert, M. Bianchi, C. Sanders, J. A. Miwa, S. Ulstrup, A. Oelsner, C. Tusche, Y.-J. Chen, D. Vasilyev, K. Medjanik, G. Brenner, S. Dziarzhytski, H. Redlin, B. Manschwetus, S. Dong, J. Hauer, L. Rettig, F. Diekmann, K. Rossnagel, J. Demsar, H.-J. Elmers, P. Hofmann, R. Ernstorfer, G. Schönhense, Y. Acremann, and W. Wurth, Time- and momentum-resolved photoemission studies using time-of-flight momentum microscopy at a free-electron laser, Rev. Sci. Instrum. 91, 013109 (2020).
  * Zhang  _et al._ [2022] H. Zhang, T. Pincelli, C. Jozwiak, T. Kondo, R. Ernstorfer, T. Sato, and S. Zhou, Angle-resolved photoemission spectroscopy, Nat. Rev. Dis. Primers 2, 54 (2022).
  * Baumgärtner  _et al._ [2022] K. Baumgärtner, M. Reuner, C. Metzger, D. Kutnyakhov, M. Heber, F. Pressacco, C.-H. Min, T. R. F. Peixoto, M. Reiser, C. Kim, W. Lu, R. Shayduk, M. Izquierdo, G. Brenner, F. Roth, A. Schöll, S. Molodtsov, W. Wurth, F. Reinert, A. Madsen, D. Popova-Gorelova, and M. Scholz, Ultrafast orbital tomography of a pentacene film using time-resolved momentum microscopy at a FEL, Nat. Commun. 13, 2741 (2022).
  * Lüftner  _et al._ [2014] D. Lüftner, T. Ules, E. M. Reinisch, G. Koller, S. Soubatch, F. S. Tautz, M. G. Ramsey, and P. Puschnig, Imaging the wave functions of adsorbed molecules, Proc. Nat. Acad. Sci. 111, 605 (2014).
  * Offenbacher  _et al._ [2015] H. Offenbacher, D. Lüftner, T. Ules, E. M. Reinisch, G. Koller, P. Puschnig, and M. G. Ramsey, Orbital tomography: Molecular band maps, momentum maps and the imaging of real space orbitals of adsorbed molecules, J. Electron Spectrosc. Relat. Phenom. 204, 92 (2015).
  * Kliuiev  _et al._ [2016] P. Kliuiev, T. Latychevskaia, J. Osterwalder, M. Hengsberger, and L. Castiglioni, Application of iterative phase-retrieval algorithms to ARPES orbital tomography, New J. Phys. 18, 093041 (2016).
  * Kliuiev  _et al._ [2018] P. Kliuiev, T. Latychevskaia, G. Zamborlini, M. Jugovac, C. Metzger, M. Grimm, A. Schöll, J. Osterwalder, M. Hengsberger, and L. Castiglioni, Algorithms and image formation in orbital tomography, Phys. Rev. B 98, 085426 (2018).
  * Jansen  _et al._ [2020] G. S. M. Jansen, M. Keunecke, M. Düvel, C. Möller, D. Schmitt, W. Bennecke, F. J. S. Kappert, D. Steil, D. R. Luke, S. Steil, and S. Mathias, Efficient orbital imaging based on ultrafast momentum microscopy and sparsity-driven phase retrieval, New J. Phys. 22, 063012 (2020).
  * Gerchberg and Saxton [1972] R. W. Gerchberg and W. O. Saxton, A Practical Algorithm for the Determination of Phase from Image and Diffraction Plane Pictures, Optik 35, 237 (1972).
  * Fienup [1982] J. R. Fienup, Phase retrieval algorithms: a comparison, Appl. Opt. 21, 2758 (1982).
  * Marchesini  _et al._ [2003] S. Marchesini, H. He, H. N. Chapman, S. P. Hau-Riege, A. Noy, M. R. Howells, U. Weierstall, and J. C. H. Spence, X-ray image reconstruction from a diffraction pattern alone, Phys. Rev. B 68, 140101 (2003).
  * Weiß  _et al._ [2015] S. Weiß, D. Lüftner, T. Ules, E. M. Reinisch, H. Kaser, A. Gottwald, M. Richter, S. Soubatch, G. Koller, M. G. Ramsey, F. S. Tautz, and P. Puschnig, Exploring three-dimensional orbital imaging with energy-dependent photoemission tomography, Nat. Commun. 6, 8287 (2015).
  * Graus  _et al._ [2019] M. Graus, C. Metzger, M. Grimm, P. Nigge, V. Feyer, A. Schöll, and F. Reinert, Three-dimensional tomographic imaging of molecular orbitals by photoelectron momentum microscopy, Eur. Phys. J. B 92, 80 (2019).
  * Kliuiev  _et al._ [2019] P. Kliuiev, G. Zamborlini, M. Jugovac, Y. Gurdal, K. v. Arx, K. Waltar, S. Schnidrig, R. Alberto, M. Iannuzzi, V. Feyer, M. Hengsberger, J. Osterwalder, and L. Castiglioni, Combined orbital tomography study of multi-configurational molecular adsorbate systems, Nat. Commun. 10, 5255 (2019).
  * Candès  _et al._ [2013] E. J. Candès, T. Strohmer, and V. Voroninski, PhaseLift: Exact and Stable Signal Recovery from Magnitude Measurements via Convex Programming, Comm. Pure Appl. Math 66, 1241 (2013).
  * Candès  _et al._ [2015] E. J. Candès, Y. C. Eldar, T. Strohmer, and V. Voroninski, Phase Retrieval via Matrix Completion, SIAM Review 57, 225 (2015).
  * [22] The use of the PhaseLift algorithm in POT was suggested to the authors by Manabu Hagiwara and Kazushi Mimura and formulation described in this work is similar to that in forthcoming paper [46].
  * Branch  _et al._ [1999] M. A. Branch, T. F. Coleman, and Y. Li, A Subspace, Interior, and Conjugate Gradient Method for Large-Scale Bound-Constrained Minimization Problems, SIAM Journal on Scientific Computing 21, 1 (1999).
  * Wales and Doye [1997] D. J. Wales and J. P. K. Doye, Global Optimization by Basin-Hopping and the Lowest Energy Structures of Lennard-Jones Clusters Containing up to 110 Atoms, The Journal of Physical Chemistry A 101, 5111 (1997).
  * Szabo and Ostlund [1996] A. Szabo and N. S. Ostlund,  _Modern Quantum Chemistry_ (Dover Publications, Inc., New York, 1996).
  * Löwdin [1950] P.-O. Löwdin, On the Non-Orthogonality Problem Connected with the Use of Atomic Wave Functions in the Theory of Molecules and Crystals, J. Chem. Phys. 18, 365 (1950).
  * Wang  _et al._ [2004] Y. L. Wang, W. Ji, D. X. Shi, S. X. Du, C. Seidel, Y. G. Ma, H.-J. Gao, L. F. Chi, and H. Fuchs, Structural evolution of pentacene on a Ag(110) surface, Phys. Rev. B 69, 075408 (2004).
  * Willenbockel  _et al._ [2013] M. Willenbockel, B. Stadtmüller, K. Schönauer, F. C. Bocquet, D. Lüftner, E. M. Reinisch, T. Ules, G. Koller, C. Kumpf, S. Soubatch, P. Puschnig, M. G. Ramsey, and F. S. Tautz, Energy offsets within a molecular monolayer: the influence of the molecular environment, New J. Phys. 15, 033017 (2013).
  * Abbasi and Scholz [2009] A. Abbasi and R. Scholz, Ab Initio Calculation of the Dispersion Interaction between a Polyaromatic Molecule and a Noble Metal Substrate: PTCDA on Ag(110), J. Phys. Chem. C 113, 19897 (2009).
  * Wießner  _et al._ [2013] M. Wießner, J. Ziroff, F. Forster, M. Arita, K. Shimada, P. Puschnig, A. Schöll, and F. Reinert, Substrate-mediated band-dispersion of adsorbate molecular states, Nat. Commun. 4, 1514 (2013).
  * Ziroff  _et al._ [2010] J. Ziroff, F. Forster, A. Schöll, P. Puschnig, and F. Reinert, Hybridization of Organic Molecular Orbitals with Substrate States at Interfaces: PTCDA on Silver, Phys. Rev. Lett. 104, 233004 (2010).
  * Metzger [2021] C. T. P. Metzger,  _Development of photoemission spectroscopy techniques for the determination of the electronic and geometric structure of organic adsorbates_ , doctoral thesis, Universität Würzburg (2021).
  * Becke [1993] A. D. Becke, A new mixing of Hartree–Fock and local density‐functional theories, J. Chem. Phys. 98, 1372 (1993).
  * Frisch  _et al._ [2016] M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson, H. Nakatsuji, X. Li, M. Caricato, A. V. Marenich, J. Bloino, B. G. Janesko, R. Gomperts, B. Mennucci, H. P. Hratchian, J. V. Ortiz, A. F. Izmaylov, J. L. Sonnenberg, D. Williams-Young, F. Ding, F. Lipparini, F. Egidi, J. Goings, B. Peng, A. Petrone, T. Henderson, D. Ranasinghe, V. G. Zakrzewski, J. Gao, N. Rega, G. Zheng, W. Liang, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, K. Throssell, J. A. Montgomery, Jr., J. E. Peralta, F. Ogliaro, M. J. Bearpark, J. J. Heyd, E. N. Brothers, K. N. Kudin, V. N. Staroverov, T. A. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. P. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, J. M. Millam, M. Klene, C. Adamo, R. Cammi, J. W. Ochterski, R. L. Martin, K. Morokuma, O. Farkas, J. B. Foresman, and D. J. Fox, Gaussian˜16 Revision C.02 (2016), gaussian Inc. Wallingford CT.
  * Kresse and Furthmüller [1996] G. Kresse and J. Furthmüller, Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set, Phys. Rev. B 54, 11169 (1996).
  * Kresse and Joubert [1999] G. Kresse and D. Joubert, From ultrasoft pseudopotentials to the projector augmented-wave method, Phys. Rev. B 59, 1758 (1999).
  * Perdew  _et al._ [1996] J. P. Perdew, K. Burke, and M. Ernzerhof, Generalized Gradient Approximation Made Simple, Phys. Rev. Lett. 77, 3865 (1996).
  * Grimme  _et al._ [2010] S. Grimme, J. Antony, S. Ehrlich, and H. Krieg, A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu, The Journal of Chemical Physics 132, 154104 (2010).
  * Methfessel and Paxton [1989] M. Methfessel and A. T. Paxton, High-precision sampling for Brillouin-zone integration in metals, Phys. Rev. B 40, 3616 (1989).
  * Monkhorst and Pack [1976] H. J. Monkhorst and J. D. Pack, Special points for Brillouin-zone integrations, Phys. Rev. B 13, 5188 (1976).
  * O’Donoghue  _et al._ [2016] B. O’Donoghue, E. Chu, N. Parikh, and S. Boyd, Conic Optimization via Operator Splitting and Homogeneous Self-Dual Embedding, Journal of Optimization Theory and Applications 169, 1042 (2016).
  * Diamond and Boyd [2016] S. Diamond and S. Boyd, CVXPY: A Python-embedded modeling language for convex optimization, Journal of Machine Learning Research 17, 1 (2016).
  * Agrawal  _et al._ [2018] A. Agrawal, R. Verschueren, S. Diamond, and S. Boyd, A rewriting system for convex optimization problems, Journal of Control and Decision 5, 42 (2018).
  * Momma and Izumi [2011] K. Momma and F. Izumi, VESTA 3 for three-dimensional visualization of crystal, volumetric and morphology data, J. Appl. Crystallogr. 44, 1272 (2011).
  * Schiefer  _et al._ [2007] S. Schiefer, M. Huth, A. Dobrinevski, and B. Nickel, Determination of the Crystal Structure of Substrate-Induced Pentacene Polymorphs in Fiber Structured Thin Films, Journal of the American Chemical Society 129, 10316 (2007).
  * Niki  _et al._ [2023] K. Niki, R. Asano, R. Sakanoue, M. Hagiwara, and K. Mimura, Photoemission orbital tomography using robust sparse phaselift (2023), arXiv:2307.12500 [cond-mat.mtrl-sci] .