Photonic band structure design using persistent homology

There is growing interest in the design of synthetic gauge fields and topological band structures for light,^1–6 motivated both by their exotic behavior and potential applications as disorder-robust photonic devices.^7,8 The most common approach has been to emulate well-known topological phases from condensed matter physics.⁹ However, recent advances in synthetic dimensions in photonic systems provide platforms for observing, for the first time, topological effects in higher dimensions.^10–12 Given the high degree of control and flexibility we now have over the design of photonic band structures and synthetic gauge fields, with many different degrees of freedom, a natural question that arises is whether there exist interesting and useful topological phenomena falling outside the well-established paradigm of topological band theory, which is concerned with the “shape” of Bloch wave eigenstates, i.e., whether they exhibit vortices or gauge discontinuities within the Brillouin zone, independent of the energy dispersion.

Controlling the shape of the energy dispersion landscape can be equally important. For example, the shape of a medium’s isofrequency contours or surfaces dictates the far-field radiation profile of localized emitters,¹³ and flat dispersion relations are highly desirable for realizing strong light–matter interactions.¹⁴ Both of these examples are based on the band structure at a specific energy. However, perfectly flat Bloch bands are idealizations never achievable in practice, and real systems inevitably exhibit losses, meaning there will be a finite energy resolution. These considerations motivate the use of “fuzzy” tools for characterizing the topology and abstract “shapes” of photonic media.

The characterization of topological properties in the presence of finite resolution or noise can be carried out using persistent homology, which is a tool from the field of topological data analysis (TDA).^15–17 Persistent homology characterizes the “shape” and connectivity of high-dimensional datasets by studying how their topological properties change as a characteristic distance scale is varied. Several studies have started to apply persistent homology and other techniques from TDA to physics problems.^18–23 For example, persistent homology has recently been used to characterize entanglement of multiqubit states,¹⁹ identify different dynamical regimes of the Bose–Hubbard model,²⁰ and recognize topological phase transitions in models of interacting spins.^21–23 To apply persistent homology, the important first step is to choose a suitable distance metric to characterize the data.

Our aim in this study is to introduce persistent homology to photonics and show how it may be useful for optimization of synthetic gauge fields and discovery of novel classes of photonic band structures, which merit more detailed studies. As a simple example, we consider the characterization of the Haldane model, a tight binding model describing a honeycomb array of waveguides or resonators in the presence of a synthetic gauge field.^24–27 We show how to use persistent homology to characterize the shape of its isoenergy contours and the connectivity of the Bloch function eigenstates in the vicinity of this energy, described by the quantum distance.^28–34 As an application of this approach, we will show how these two features can be used to optimize the lifetimes and radiation profiles of embedded quantum emitters.

The outline of this article is as follows: In Sec. II, we provide a brief introduction to the technique of persistent homology and discuss how it may be applied to band structure optimization by using suitable distance metrics. Section III applies persistent homology to analyze the band structure of the Haldane model, identifying parameter ranges exhibiting “moat band”^35,36 and multivalley dispersion relations.^37,38 In Sec. IV, we show how the considered distance metrics can be used to predict properties of quantum emitters embedded in the lattice. Section V concludes with discussion of possible future directions and applications of persistent homology in photonics.

The following discussion is quite condensed and not intended to be a comprehensive introduction to persistent homology. For a more pedagogical and complete introduction, we recommend the excellent Ref. 18 aimed at physicists, or other recent articles applying persistent homology to characterize quantum and condensed matter systems.^19–22 

Persistent homology is a technique used to characterize sets of data $x1,x2,…,xN$ equipped with some measure of distance between them d(x_i, x_j) ≥ 0. For example, the data x_i could correspond to particle positions in real space characterized by the Euclidean metric d(x_i, x_j) = |x_i − x_j|, or image data with x_i encoding the locations and intensities of the individual pixels. Given the data, one can construct subsets known as k-simplices, where each k-simplex is a set of k + 1 of the points; 0-simplices are individual points {x_i}, 1-simplices are edges connecting pairs of points {x_i, x_j}, 2-simplices are areas enclosed by three points {x_i, x_j, x_k}, and so on, as illustrated in Fig. 1.

To characterize the shape of the data, we construct sets of simplices using the distance measure d to form simplicial complexes. A k-complex is a closed set of simplices of dimension up to k. By closed, we mean that, if a simplex belongs to the complex, then its edges are also in the complex. For example, if an edge {x_i, x_j} belongs to the complex, then its constituent points {x_i}, {x_j} also belong to the complex. In the case of Euclidean data, a common choice is the Vietoris–Rips complex $Vε$⁠, which includes all k-simplices whose constituent points all have a pairwise distance less than ε from each other.

By counting the number of simplices of each dimension in the complex, we can determine its abstract “shape.” For example, a complex consisting of three points and two edges describes a line, while a complex with three points and three edges forms a loop. Formally, the existence of topologically nontrivial structures of dimension k in the simplicial complex is determined by computing its kth Betti number B_k, which corresponds to the number of k-dimensional “holes” in the complex.

Adding a single vertex or edge to the complex can drastically alter its topological properties encoded by the Betti numbers. The key innovation provided by persistent homology, which makes it useful for characterizing real world systems with noise, is to study how the topology of the simplicial complex varies with some characteristic scale ε, in effect studying the shape of the data over a range of scales, known as filtration. Features that only persist for a small range of scales can be attributed to noise and discarded, while those that persist for a large range of scales are robust and meaningful.

There are several possible choices for the filtration, depending on the data to be analyzed. For the case of point cloud data, one can use the Vietoris–Rips complex $Vε$⁠. However, this becomes extremely time-consuming to compute for large numbers of data points. For point cloud data, there are other, more efficient complexes such as the α complex.^15,16 For image data (pixels on a regular grid), another choice is the sublevel set filtration constructed out of all neighboring pixels with intensities less than (or greater than) some threshold value ε, which characterizes the number and shape of intensity minima (or maxima).

Given a dataset, distance metric, and choice of filtration, the persistent homology of the data can be computed using existing software libraries.^39–41 As long as the dataset is not too large, low dimensional topological features can be computed relatively quickly. However, the computation time grows rapidly for high-dimensional features, generally requiring pre-processing to reduce the dataset size.

Let us now return to physics. In order to apply persistent homology to physical systems, we first need to identify the relevant dataset {x_i} and a suitable distance metric to characterize the data. Several different approaches have recently been proposed in the physics literature.^18–22 Here, our aim is to use persistent homology to characterize energy bands of periodic photonic media.

We will consider, for simplicity, Bloch wave spectra of photonic waveguide lattices obtained from the Schrödinger eigenvalue problem,⁷ 

where $Ĥ(k)$ is the lattice’s Bloch Hamiltonian, k is the Bloch momentum (restricted to the first Brillouin zone), |u_n(k)⟩ is the Bloch function, with the bra-ket notation encoding the internal (e.g., sublattice or spin) degrees of freedom, n is the band index, and ω_n is the frequency of the mode. We note that Maxwell’s equations for electromagnetic modes of periodic photonic media can be cast into a similar form to Eq. (1), making our approach equally applicable to the design and optimization of photonic crystals.

Numerical solutions of Eq. (1) typically sample k on a regular grid in the Brillouin zone, providing a discrete set of data (the eigenvectors and eigenvalues) to which we can apply persistent homology. Keeping with the spirit of successful applications of machine learning to physics, our aim is to employ simple distance metrics with well-established importance that can be easily interpreted.

First, we consider classifying the energy eigenvalues of a band ω(k) using a similar approach to that used for characterizing image data. That is, we sample ω on a k mesh spanning the Brillouin zone. This mesh can be a regular grid respecting the lattice symmetries (e.g., triangular, cubic, etc.), or even a random grid, provided a sufficiently large number of points are sampled. We then define a filter function f(k) = |ω(k) − ω₀|, where ω₀ is some reference frequency of interest. Persistent homology can characterize the shape of the set of momenta k, whose Bloch mode frequencies are within ε_ω of the reference frequency ω₀, i.e.,

In particular, persistent homology can tell us the range of frequency resolutions ε_ω over which topological features of the lattice’s isoenergy contours persist.

An important consideration is that we inevitably sample k on a discrete grid of length N; hence, there will be noise induced by the discretization. To be precise, neighboring grid points will have an energy difference δω ≈ $v$_g · δk, where $v$_g = ∇_kω is the local wave group velocity and δk = 2π/N is the grid spacing (using units such that the lattice period a = 1). Thus, topological features of the isoenergy surfaces are only meaningful if they persist over a range of energies larger than δω = max_κ$v$_g/δk. As a crude estimate, if Δ is the width of the band of interest, the group velocity will have a scale $vg∼Δ2π$⁠. Therefore, as a rule of thumb, one should look for features with a minimum persistence of Δ/N to avoid discretization errors; in the following, we will use δω = 2max_κv_g/N as the persistence cutoff.

The calculation proceeds as follows: (1) Choose a reference energy ω₀ and grid size N. (2) Compute the persistent homology of the isoenergy contours, using the energy resolution ϵ_ω as the filter function. (3) Discard features with persistence less than δω. (4) Sum the number of remaining features existing at each ϵ_ω of interest to obtain the Betti numbers, which characterize the shape of the isoenergy surface at that resolution, i.e., its connectivity and number of holes.

As an example of this procedure, we consider the energy band of a simple square lattice shown in Fig. 2(a), described by the dispersion relation

where J₁ and J₂ are the nearest and next-nearest neighbor hopping strengths, respectively. We take the reference energy to be ω₀ = 0, close to the band center, and J₂ = 0.2J₁. The dispersion ω(k) is shown in Fig. 2(b).

Applying persistent homology to this two-dimensional dispersion “image,” we obtain the energy resolution-dependent Betti numbers B_0,1(ϵ_ω) shown in Fig. 2(c). For small ε, there is a sharp peak in the number of distinct clusters B₀, which is induced by discreteness of the momentum grid and highly sensitive to the chosen grid size N. From this plot, we see that thresholding, i.e., consideration of features only with a minimum persistence, is essential to obtaining meaningful Betti numbers.

Applying the minimum persistence threshold in Fig. 2(d), we obtain two “phases” with different first Betti numbers B₁, in which we count the number of loops (or holes) in the data. Figure 2(e) illustrates the first phase, considering modes with energies |ω| < 2J₁, where 2J₁ is the chosen energy resolution. B₀ = 1 indicates the modes form a single connected cluster, as can be readily verified. B₁ = 3 means that the image hosts three inequivalent non-contractible loops, as shown. We note that one needs to be careful with counting non-contractible loops here, due to the periodic boundary conditions of the Brillouin zone.

As the energy cutoff ϵ_ω is increased, the band minimum at ω ≈ −3J₁ becomes included in the filtration, destroying one of the holes and resulting in the image shown in Fig. 2(f). The loss of a hole corresponds to B₁ decreasing by 1. Note that the removal of the remaining hole in the center of the image does not affect B₁; due to the Brillouin zone’s periodic boundary conditions, any loop encircling this hole can be continuously deformed to a point by first expanding it to the Brillouin zone edges at k_x = ±π to create a pair of loops and then annihilating the pair at the k_y = ±π Brillouin zone edge.

Persistent homology allows us to also classify more complicated objects than energy band structure “images.” In particular, topological phases of light and matter are associated with non-trivial topological properties of the Bloch function eigenstates |u_n(k)⟩, which are usually computed using the Berry connection or curvature.⁷ Here, we will consider the abstract shape of the Bloch functions characterized by the quantum distance d,

The quantum distance between Bloch functions at infinitesimally separated momenta corresponds to the quantum metric tensor $ĝ$⁠, with components

The study of interesting observables related to properties of the quantum distance and quantum metric currently attracts broad interest.^28–34 For example, in graphene, the large quantum distance between counter-propagating modes at the Dirac cones is responsible for the suppression of backscattering in the presence of long range (valley-conserving) disorder;²⁸ a small quantum distance means there can be efficient scattering between the two states by scalar (e.g., spin-preserving) disorder.

We use the quantum distance to define a graph, with each vertex corresponding to a Bloch function with some momentum k. To apply persistent homology, we define two points k and k′ to be connected when d(k, k′) < ϵ_d. The Betti numbers of the resulting graph can then be used to characterize the “shape” of the eigenstates, such as the existence of loops or multiple disconnected clusters. One can either consider how the number of features changes as a function of ϵ_d, or consider the features at a particular ϵ_d.

As with the energy eigenvalues, sampling of the Brillouin zone on a discrete grid will induce discretization noise for small quantum distance cutoffs ϵ_d. For neighboring grid points, $d(k,k+δk)≈δk⋅ĝ(k)⋅δk$⁠. Therefore, to avoid discretization errors, one should consider features that persist for $≳2πNĝ$⁠.

We note that under this definition, the existence of non-trivial shapes of the energy eigenvalues is a necessary (but not sufficient) condition for the Bloch functions to have non-trivial shapes, since the Bloch functions are typically continuous functions of k. Furthermore, while the energy eigenvalues share the periodicity of the reciprocal lattice, i.e., ω_n(k) = ω_n(k + G) for any reciprocal lattice vector G, the quantum distance generally does not share this periodicity. Thus, d(k, k + G) can be nonzero. In the following, we will compute topological features of the Bloch functions in the first Brillouin zone.

As an example, we consider the connectivity of the Bloch function eigenstates in the one-dimensional Su–Schrieffer–Heeger (SSH) model shown in Fig. 3, which is described by the Bloch Hamiltonian,⁷ 

where J_1,2 are coupling strengths and $σ^i$ are Pauli matrices. The eigenvectors of Eq. (6) are

and lie on the unit circle. The eigenvectors span the entire circle in the nontrivial phase J₂ > J₁, and are limited to −π/2 < θ < π/2 inthe trivial phase. The quantum distance between eigenstates at k and k′ can be evaluated exactly as

attaining its maximal value of 1 only in the nontrivial phase. We sample the Brillouin zone on a uniform grid k = k_n = 2πn/N, using N = 30 grid points. Figures 3(c) and 3(d) show the Betti numbers obtained from the two phases, without using any persistence threshold. Similar to the case of the energy eigenvalue images, discretization-induced noise appears for small filtration values ϵ_d as N-dependent peaks in B₀.

To systematically remove this noise, we introduce a persistence threshold based on the quantum metric, which in the SSH model takes the form

Using a threshold of $2.5∑n2πN2gkk(kn)$ (i.e., proportional to the average quantum distance), we obtain smooth curves distinguishing the trivial and nontrivial phases in Figs. 3(e) and 3(f). In particular, the winding of the Bloch functions in the nontrivial phase is characterized by B₁ = 1, which persists over a large range of scales ϵ_d. We note that as ϵ_d → 1, each point becomes connected to all others, destroying the hole and resulting in B₁ = 0. Thus, to characterize the shape of the Bloch functions, one should use an intermediate value of ϵ_d.

We now apply the above measures to the classification of the Haldane model shown in Fig. 4(a),^24–27 described by the Bloch wave Hamiltonian

where $δ1=(0,1/3)$⁠, $δ2=12(−1,−1/3)$⁠, and $δ3=12(1,−1/3)$ are displacements between neighboring unit cells, a₁ = (1, 0), $a2=12(1,3)$⁠, and $a3=12(−1,3)$ are the lattice vectors, J₁, J₂ ≥ 0 are nearest- and next-nearest-neighbour hopping strengths, M is the sublattice detuning, and ϕ parameterizes the staggered next-nearest-neighbour hopping phases.

A previous study showed that when M = ϕ = 0, the Haldane model can exhibit “moat bands,” i.e., degenerate band minima forming closed lines in the Brillouin zone.³⁵ The existence condition for moat bands when M = ϕ = 0 was obtained by rewriting $Ĥ(k)$ as

which implies its energy eigenvalues take the form $ω1,2=±J1Tx2+Ty2+J2(Tx2+Ty2)$⁠. Consequently, the eigenvalues exhibit a degenerate minimum along the contour $Tx2+Ty2=J1/(2J2)$⁠. Since $Tx2+Ty2≤3$⁠, moat bands can only exist when the next-nearest neighbor coupling strength is sufficiently strong, i.e., J₂ > J₁/6. When a moat band is formed, the ground state of hard-core bosons at low densities spontaneously breaks time-reversal symmetry.³⁵ 

When M, ϕ ≠ 0, the simple condition for the existence of a moat band no longer holds. We apply the persistent homology approach discussed in Sec. II B to systematically determine the fate of the moat band for M, ϕ ≠ 0. To do so, we set the reference energy as ω₀ = min_k[ω₁(k)], considering the shape of the lowest energy eigenvalues of $Ĥ(k)$⁠. We sample the Brillouin zone on a grid of N × N = 30 × 30 points and consider features at a fixed frequency resolution ϵ_ω = 0.1J₁. Figures 4(b)–4(f) show examples of the band structure in different phases identified using persistent homology, exhibiting different shapes of the lowest energy modes characterized by the Betti numbers plotted in Fig. 5.

First, when M = ϕ = 0, we observe that the persistent homology approach accurately reproduces the features of the energy minima obtained in Ref. 35. Note that the threshold for the emergence of a single loop is slightly larger than the analytical threshold J₁ = J₁/6, due to the finite frequency resolution we consider, which results in a filled disk for small loop radii. Figure 4(b) shows an example of a single-looped energy minimum (B₁ = 1), which persists for nonzero ϕ or M. As J₂ is increased, the loop expands, generating additional loops once it intersects the edges of the Brillouin zone, characterized by Betti numbers B₁ = 4 [Fig. 4(c)] or 2 [Fig. 4(e)], depending on the parameter values. For large J₂, the single connected contour splits into two isolated ring-shaped energy minima centered at the two K points. The differing behavior under increasing ϕ and M is noteworthy; increasing ϕ always results in the destruction of the loops and their replacement with a pair of isolated, point-like energy minima [Fig. 4(d)], whereas the loops can persist in the presence of moderate M.

Interestingly, by fine-tuning a pair of parameters (J₂, ϕ) or (J₂, M), one can even obtain a “phase transition point” at which multiple Betti numbers B_0,1 intersect [Fig. 4(f)]. At this point, it is possible to obtain a variety of shapes of the low energy modes by applying small perturbations.

It is also noteworthy that, for the chosen energy resolution ϵ_ω = 0.1J₁, the nontrivial low energy mode shapes emerge when the lower band is much flatter than the upper band, making the shapes difficult to directly identify from the images of the full band structure in Fig. 4. One could alternatively normalize the energy resolution by the width of the band of interest, or consider features at other energies. For example, if ω₀ is set to the band gap, the nonzero overlap between the upper and lower bands in Figs. 4(c), 4(d), and 4(f) can generate multiple disconnected loops with B_0,1 > 1.

Next, in Fig. 6, we study the connectivity of the low energy Bloch functions characterized by the quantum distance. We consider features at ϵ_d = 0.5, motivated by the SSH model example discussed in Sec. II C. The Bloch function connectivity provides additional information about the shape of the band properties, complementary to the shape of the low energy eigenvalues in Fig. 5. Since the Bloch functions are continuous functions of the momentum k, they can only form disconnected clusters if the isoenergy contours also form disconnected clusters. This explains why the B₀ = 2 states appearing in Fig. 6(a) are fully contained in the B₀ = 2 region of Fig. 5(a). Similarly, for the case of B₁, isoenergy loops are necessary for loops of the Bloch functions to appear.

Moderate ϕ and M have qualitatively different effects on the shape of the Bloch functions; moderate ϕ results in the creation of two disconnected clusters of eigenvectors, while under increasing M, the Bloch functions only form a single cluster. To understand these differing behaviors, we note that at large J₂, the band minima are centered at the K points K_± = (±4π/3, 0). A long wavelength expansion of Eq. (11) yields the continuum Hamiltonian

where M and ϕ have qualitatively different effects on the Dirac mass $M∓J23sin⁡ϕ$⁠. When ϕ = 0, the two valleys have the same mass, meaning that there is a small quantum distance between them; hence, the low energy Bloch functions form a single connected cluster. On the other hand, ϕ creates two valleys with opposite masses, and hence eigenvectors reside on opposite sublattices with a large quantum distance between them, corresponding to two clusters.

Even though the correspondence between the Betti numbers and the more familiar band topological invariants such as the Chern number C is not exact, the latter may be used to guide the design of Bloch function clusters and loops. For example, in trivial C = 0 phases, the Bloch functions of a band can be continuously deformed to a point such that the quantum distance vanishes.⁷ Hence, the Bloch functions are likely to form a single cluster at any energy resolution. On the other hand, such a deformation is impossible in non-trivial C ≠ 0 phases and there always exists Bloch functions separated by a large quantum distance, which enables the formation of distinct clusters. In the case of the Haldane model, phases with non-trivial Chern number C occur when the two valleys have opposite signs of the Dirac mass.

In this section, we will show how our chosen metrics can be used to predict the behavior of interesting physical observables. Namely, we will demonstrate that the shape and connectivity of the energy band minima affect the lifetime of localized quantum emitters coupled to the lattice and resonant with the band edge, by imposing constraints on the effective emitter-Bloch wave coupling strength.

We consider the coupling of N_e localized quantum emitters to the Haldane lattice, described by the Fano–Anderson model (ℏ = 1)^42–45 

where $ĉα$ annihilates a bound state with energy ω_α and non-radiative decay rate γ_α, $ûn(k)$ annihilates a Bloch wave in band n with momentum k, and g_α,n(k) is the emitter-Bloch wave coupling strength,

which we have decomposed into an energy scale J_e and a dimensionless part $|⟨un(k)|gα(k)⟩|=1−dα,n(k)$⁠, where d_α,n(k) is the quantum distance between the emitter and Bloch wave polarizations at momentum k. For simplicity, we focus on the single excitation subspace, in which the emitter and Bloch wave probability amplitudes evolve according to

We assume that only the emitters are excited at t = 0, i.e., u_n(k, 0) = 0, and applying the weak coupling and Markov approximations (see Refs. 42 and 43), we arrive at an effective evolution equation for the emitter amplitudes,

where $Hα,β=(ωα−iγα)δα,β+Δα,β$⁠,

describes the emitters’ coupling and self-energy mediated by the lattice, and we have dropped the summation over the band index n, as only the band energetically closest to the emitters will be relevant in the weak coupling regime. Note that, under the weak coupling approximation, J_e only sets the scale of Δ_α,β and is otherwise irrelevant.

When the emitter frequencies are close to the lower band edge, the integral in Eq. (21) will be dominated by modes with energies within the emitter linewidth |ω₁(k) − ω_β| < γ_β such that Δ_α,β becomes sensitive to the shape and connectivity of the eigenvalues and Bloch functions within this energy resolution. For example, when the emitter couplings g_α,1 ≈ 0 at the high symmetry points of the Brillouin zone, the emitter-lattice coupling will be strongly suppressed unless the band minimum forms a “moat” band.

We illustrate this in Fig. 7 for the case of a single emitter tuned to the band minimum with non-radiative decay rate γ_e = 0.1J₁ and coupled to a six-site plaquette to form a “giant atom,”^44,45 described by the coupling

Figures 7(a) and 7(b) show the radiative decay rate Γ = −Im[Δ_e,e] as a function of J₂, M, and ϕ. The radiative decay rate is smallest when J₂ = 0, corresponding to a point-like band minimum at k = 0. The peaks in the decay rate coincide reasonably well with the formation of the moat bands. To show that this enhancement of Γ is sensitive to the shape of the low energy modes and not merely determined by the number of resonant modes, we also show in Figs. 7(c) and 7(d) the normalized decay rate $Γ/N$⁠, where

is the |g₁(k)⟩-independent part of Eq. (21). The normalized decay rate is largest when there are multiple loops, i.e., B₁ = 2 or 4, which optimizes the trade-off between emitter-Bloch wave mode detunings ω₁(k) − ω_e and the emitter-Bloch wave coupling ⟨u_n(k)|g_e(k)⟩.

We have shown how a method from the field of topological data analysis, namely, persistent homology, can be used to identify novel classes of periodic lattices by systematically classifying the abstract shapes of their energy dispersion and Bloch functions. By using physically motivated choices for the minimal persistence, we were able to eliminate discretization-induced noise in the topological features, thereby enabling reliable high-throughput classification of photonic band structures. The approach may be used to save many graduate student-hours, for example, if it is necessary to compute and check thousands of band structures to find an optimal lattice design subject to certain constraints. Moreover, the method is highly flexible and can be readily generalized to higher dimensional settings,^10,11 complex lattices with large unit cells or overlapping bands such as Moiré superlattices,^46,47 continuous media such as photonic crystals,⁴⁸ and nonlinear or strongly interacting topological wave media.⁸ 

For the simple tight binding models we considered, the small parameter space could be covered through brute-force scanning. Persistent homology is compatible with gradient-based optimization approaches, so an interesting direction for future research will be to use gradient-based methods to find specific energy band topologies of interest given constraints in the lattice design.

In contrast to some other recent studies employing persistent homology in physics,^18–23 rather than focusing on the construction of persistence diagrams and their similarity measures, we focused on distance metrics directly related to physical observables of interest, such as the energy resolution. This allowed us to reliably separate robust features of the band structure from the topological noise induced by numerical discretization. Future work may optimize the behavior of other types of observables using persistent homology, study connections between the Betti numbers, and the usual invariants of topological band theory such as the Chern number and explore applications of persistent homology to other areas of photonics such as imaging. We also hope this approach will also attract interest in branches of physics beyond photonics, such as the design and classification of novel condensed matter and quantum systems.

The data that support the findings of this study are available from the corresponding author upon reasonable request.

This research was supported by the National Research Foundation, Prime Ministers Office, Singapore, the Ministry of Education, Singapore, under the Research Centres of Excellence program, and the Polisimulator project co-financed by Greece and the EU Regional Development Fund.