### Pattern Matching in large Graph Transformation Systems

Chapter 4

# Pattern Matching in large Graph Transformation Systems

To our knowledge, the first practical proposal for a GTS-based quantum compiler was presented in Xu, 2022Mingkuan Xu, Zikun Li, Oded Padon, Sina Lin, Jessica Pointing, Auguste Hirth, Henry Ma, Jens Palsberg, Alex Aiken, Umut A. Acar and Zhihao Jia. 2022. Quartz: Superoptimization of Quantum Circuits. In Proceedings of the 43rd ACM SIGPLAN International Conference on Programming Language Design and Implementation, June 2022. Association for Computing Machinery, 625--640. doi: 10.1145/3519939.3523433 and then refined in Xu, 2023Amanda Xu, Abtin Molavi, Lauren Pick, Swamit Tannu and Aws Albarghouthi. 2023. Synthesizing Quantum-Circuit Optimizers. Proceedings of the ACM on Programming Languages 7, PLDI (June 2023, 835--859). doi: 10.1145/3591254. In these, the set of possible graph transformations is obtained by exhaustive enumeration. Using SAT solvers and fingerprinting techniques, the set of all small programs up to a certain size can be generated ahead of time and clustered into disjoint partitions of equivalent programs. This concisely expresses every possible peephole optimisation up to the specified size: for every small enough subset of operations of an input program, its equivalence class can be determined. Any replacement of that set of operations with another program in the same equivalence class is a valid transformation and, thus, a potential peephole optimisation. Transformation systems on minIR graphs based on equivalence classes were formalised in section 3.4.

First results of this approach are promising. Xu, 2022Mingkuan Xu, Zikun Li, Oded Padon, Sina Lin, Jessica Pointing, Auguste Hirth, Henry Ma, Jens Palsberg, Alex Aiken, Umut A. Acar and Zhihao Jia. 2022. Quartz: Superoptimization of Quantum Circuits. In Proceedings of the 43rd ACM SIGPLAN International Conference on Programming Language Design and Implementation, June 2022. Association for Computing Machinery, 625--640. doi: 10.1145/3519939.3523433 demonstrated that optimisation performance improves markedly with larger sets of transformation rules. Such workloads however rely heavily on pattern matching, the computational task that identifies subgraphs on which transformation rules apply. In Xu, 2022Mingkuan Xu, Zikun Li, Oded Padon, Sina Lin, Jessica Pointing, Auguste Hirth, Henry Ma, Jens Palsberg, Alex Aiken, Umut A. Acar and Zhihao Jia. 2022. Quartz: Superoptimization of Quantum Circuits. In Proceedings of the 43rd ACM SIGPLAN International Conference on Programming Language Design and Implementation, June 2022. Association for Computing Machinery, 625--640. doi: 10.1145/3519939.3523433 and Xu, 2023Amanda Xu, Abtin Molavi, Lauren Pick, Swamit Tannu and Aws Albarghouthi. 2023. Synthesizing Quantum-Circuit Optimizers. Proceedings of the ACM on Programming Languages 7, PLDI (June 2023, 835--859). doi: 10.1145/3591254, pattern matching is carried out separately for each pattern. This becomes a significant bottleneck for large rule sets. In Xu, 2022Mingkuan Xu, Zikun Li, Oded Padon, Sina Lin, Jessica Pointing, Auguste Hirth, Henry Ma, Jens Palsberg, Alex Aiken, Umut A. Acar and Zhihao Jia. 2022. Quartz: Superoptimization of Quantum Circuits. In Proceedings of the 43rd ACM SIGPLAN International Conference on Programming Language Design and Implementation, June 2022. Association for Computing Machinery, 625--640. doi: 10.1145/3519939.3523433, performance peaks at around 50,000 transformation rules, after which the additional overhead from pattern matching becomes dominant, deteriorating the compilation results.

In this chapter, we solve these scaling difficulties by presenting an algorithm for pattern matching on minIR graphs that uses a pre-computed data structure to return all pattern matches in a single query. The set of transformation rules is directly encoded in this data structure. After a one-time cost for construction, pattern-matching queries can be answered in running time independent of the number of rules in the transformation system.

The asymptotic complexity results presented in this chapter depend on some simplifying assumptions on the properties that the pattern graphs and embeddings must satisfy. This represents a restriction on the generality of minIR graphs, but we do not find that they restrict the usefulness of the result in practice. As discussed in section 4.7, none of these assumptions are required in practice for the implementation. We have not observed any impact on performance when the imposed constraints are lifted, so we conjecture that at least some of these assumptions can be relaxed and our results generalised.

After a discussion of related work in section 4.1, Section 4.2 presents the assumptions that we are making in detail, along with some relevant definitions for the rest of the chapter. Sections 4.3, 4.4, and 4.5 present the core ideas of our approach, respectively introducing: a reduction of minIR graphs to equivalent trees, a canonical construction for the tree reduction and an efficient way to enumerate all possible subtrees of a graph. We also prove bounds on the size and number of the resulting trees.

In section 4.6, we introduce a pre-computation step and show that the pattern-matching problem reduced to tree structures can be solved using a prefix tree-like automaton that is fixed and pre-computed for a given set of patterns. Combining the automaton construction with bounds from section 4.5 leads to the final result. We conclude in section 4.7 with benchmarks on a real-world dataset of 10000 quantum circuits, obtaining a 20x speedup over a leading C++ implementation of pattern matching for quantum circuits.

## 4.1. Related work

Our proposed solution can be seen as a specialisation of RETE networks Forgy, 1982Charles L. Forgy. 1982. Rete: A fast algorithm for the many pattern/many object pattern match problem. Artificial Intelligence 19, 1 (Septempter 1982, 17--37). doi: 10.1016/0004-3702(82)90020-0 Varró, 2013Gergely Varró and Frederik Deckwerth. 2013. A Rete Network Construction Algorithm for Incremental Pattern Matching and derivatives Ian, 2003Wright Ian and James A. R. Marshall. 2003. The execution kernel of RC++: RETE*, a faster RETE with TREAT as a special case. International Journal of Intelligent Games and Simulation 2, 1 (Feb 2003, 36-48) Armstr., 2014Dylan Armstrong. 2014. Memory Efficient Stream Reasoning onResource-Limited Devices. PhD Thesis. University of Dublin, Trinity College Mirank., 1987D. P. Miranker. 1987. TREAT: a new and efficient match algorithm for AI production systems. PhD Thesis. Columbia University, 2960 Broadway New York, NYUnited States to the case of graph pattern matching. The additional structure obtained from restricting our considerations to graphs results in a simplified network design that allows us to derive worst-case asymptotic runtime and space bounds that are polynomial in the parameters relevant to our use case1 – overcoming a key limitation of RETE.

Another well-studied application of large-scale pattern matching is in the context of stochastic biomolecular simulations Sneddon, 2010Michael W Sneddon, James R Faeder and Thierry Emonet. 2010. Efficient modeling, simulation and coarse-graining of biological complexity with NFsim. Nature Methods 8, 2 (December 2010, 177--183). doi: 10.1038/nmeth.1546 Bachman, 2011John A Bachman and Peter Sorger. 2011. New approaches to modeling complex biochemistry. Nature Methods 8, 2 (January 2011, 130--131). doi: 10.1038/nmeth0211-130, particularly the Kappa project Danos, 2004Vincent Danos and Cosimo Laneve. 2004. Formal molecular biology. Theoretical Computer Science 325, 1 (Septempter 2004, 69--110). doi: 10.1016/j.tcs.2004.03.065. Stochastic simulations depend on performing many rounds of fast pattern-matching for continuous Monte Carlo simulations Yang, 2008Jin Yang, Michael I. Monine, James R. Faeder and William S. Hlavacek. 2008. Kinetic Monte Carlo method for rule-based modeling of biochemical networks. Physical Review E 78, 3 (Septempter 2008, 031910). doi: 10.1103/physreve.78.031910. However, unlike our use case, the procedure typically does not need to scale well to a large number of patterns. In Danos, 2007Vincent Danos, Jérôme Feret, Walter Fontana and Jean Krivine. 2007. Scalable Simulation of Cellular Signaling Networks, Danos et al. introduced a pre-computation step to accelerate matching by establishing relations between patterns that activate or inhibit further patterns. This idea was later expanded upon and formalised in categorical language in Boutil., 2017Pierre Boutillier, Thomas Ehrhard and Jean Krivine. 2017. Incremental Update for Graph Rewriting. The ideas presented in Boutil., 2017Pierre Boutillier, Thomas Ehrhard and Jean Krivine. 2017. Incremental Update for Graph Rewriting are similar to ours; their formalism has the advantage of being more general but does not present any asymptotic complexity bounds and suffers from similar worst-case complexities as RETE.

A similar problem has also been studied in the context of multiple-query optimisation for database queries Sellis, 1988Timos K. Sellis. 1988. Multiple-query optimization. ACM Transactions on Database Systems 13, 1 (March 1988, 23–52). doi: 10.1145/42201.42203 Ren, 2016Xuguang Ren and Junhu Wang. 2016. Multi-Query Optimization for Subgraph Isomorphism Search. Proceedings of the VLDB Endowment 10, 3 (November 2016, 121–132). doi: 10.14778/3021924.3021929, but it has limited itself to developing caching strategies and search heuristics for specific use cases. Finally, using a pre-compiled data structure for pattern matching was already proposed in Messmer, 1999Bruno T. Messmer and Horst Bunke. 1999. A decision tree approach to graph and subgraph isomorphism detection. Pattern Recognition 32, 12 (December 1999, 1979--1998). doi: 10.1016/S0031-3203(98)90142-X. However, with a n Θ(m)n^{\Theta(m)}n Θ(m) space complexity – n nn is the input size and m mm the pattern size – it does not scale to large input graphs, even for small patterns.

1. RETE networks have been shown to have exponential worst-case space (and thus time) complexity Rakib, 2018Abdur Rakib and Ijaz Uddin. 2018. An Efficient Rule-Based Distributed Reasoning Framework for Resource-bounded Systems. Mobile Networks and Applications 24, 1 (October 2018, 82--99). doi: 10.1007/s11036-018-1140-x, although performance in practical use cases can vary widely Uddin, 2016Ijaz Uddin, Hafiz Mahfooz Ul Haque, Abdur Rakib and Mohamad Rafi Segi Rahmat. 2016. Resource-Bounded Context-Aware Applications: A Survey and Early Experiment. ↩︎

## 4.2. Preliminaries and simplifying assumptions

For simplicity, we will throughout consider minIR graphs that admit a type system Σ\SigmaΣ, though most of the results can also be adapted to other graph domains. We will write T T T for the types of Σ\SigmaΣ (i.e. its values) and Γ\GammaΓ for the operation types (i.e. its edges).

#### Linear paths and operation splitting #

An operation type ν∈Γ\nu \in \Gamma ν∈Γ in the type system Σ\SigmaΣ is a hyperedge. Its endpoints

src(ν)=src 1(ν)⋯src s(ν) and tgt(ν)=tgt 1(ν)⋯tgt t(ν) \textit{src}(\nu) = \textit{src}_1(\nu) \cdots \textit{src}_s(\nu) \textrm{ and } \textit{tgt}(\nu) = \textit{tgt}_1(\nu) \cdots \textit{tgt}_t(\nu) src(ν)=src 1​(ν)⋯src s​(ν) and tgt(ν)=tgt 1​(ν)⋯tgt t​(ν)

are strings of data types that define the input and output signature of the operation ν\nu ν. We can refer to the set of all hyperedge endpoints of ν\nu ν using the string indices I d x(⋅)\mathrm{Idx}(\cdot)Idx(⋅) (⊔\sqcup⊔ denotes the disjoint set union):

P ν=I d x(src(ν))⊔I d x(tgt(ν))={1,…,s}⊔{1,…,t}.P_\nu = \mathrm{Idx}(\textit{src}(\nu)) \sqcup \mathrm{Idx}(\textit{tgt}(\nu)) = \{1, \dots, s\} \sqcup \{1, \dots, t\}.P ν​=Idx(src(ν))⊔Idx(tgt(ν))={1,…,s}⊔{1,…,t}.

Fix a partition of P ν P_\nu P ν​ into disjoint pairs

P ν={p 1, p 1′} ⊔ {p 2, p 2′} ⊔ ⋯ ,P_\nu = \{p_1, p_1'\} \,\sqcup\, \{p_2, p_2'\} \,\sqcup\, \cdots, P ν​={p 1​,p 1′​}⊔{p 2​,p 2′​}⊔⋯,

where the last set of the partition may be a singleton if ∣P ν∣|P_\nu|∣P ν​∣ is odd. For every ν∈Γ\nu \in \Gamma ν∈Γ, we then define f=⌈∣P o∣/2⌉f = \lceil |P_o| / 2 \rceil f=⌈∣P o​∣/2⌉ new split operation types ν 1,…,ν f\nu_1, \dots, \nu_f ν 1​,…,ν f​, each with two endpoints: the i ii-th operation type ν i\nu_i ν i​ has endpoints p i p_ip i​ and p i′p_i'p i′​ in P ν P_\nu P ν​. For every operation o∈E o \in Eo∈E of type ν\nu ν, we can then split o oo into f f f operations o 1,…,o f, o_1, \dots, o_f,o 1​,…,o f​, each of arity 1 or 2 and of types ν 1,…,ν f\nu_1, \dots, \nu_f ν 1​,…,ν f​ respectively. We will refer to the graph transformation that replaces an operation o oo in a minIR graph with the operations o i o_io i​ for 1⩽i⩽f 1 \leqslant i \leqslant f1⩽i⩽f as operation splitting.

It is important to note that the splitting of an operation o oo is unique and given by the type of o oo, and thus invariant under (typed) morphisms: there is a morphism P→G P \to G P→G of a pattern into a graph if and only if there is a morphism P′→G′P' \to G'P′→G′ from the split pattern P′P'P′ into the split graph G′G'G′.

A transformation rule splitting an operation with 3 sources and 2 targets. The choice of endpoint partition made here, obtained by pairing the i ii-th use with the i ii-th define, is arbitrary but convenient for quantum gates as they correspond to the input and output values of a same qubit.

The endpoint partitions P ν P_\nu P ν​ also define linear paths. Two values v, v′v, v'v, v′ in a minIR graph are on the same linear path if there are values u 1,…,u k u_1, \dots, u_ku 1​,…,u k​ with v=u 1 v = u_1 v=u 1​ and v′=u k v' = u_k v′=u k​ such that u i u_iu i​ is connected to u i+1 u_{i+1}u i+1​ through an operation o oo and they correspond to the same pair of endpoints in the endpoints partition (i.e. the indices of P t y p e(o)P_{type(o)}P t y pe(o)​ correspond to values u i u_iu i​ and u i+1 u_{i+1}u i+1​ in src(o)⊔tgt(o)\textit{src}(o) \sqcup \textit{tgt}(o)src(o)⊔tgt(o)).

#### Linearity assumption and rigidity #

Recall that in Definition 3.2, V src V_\textit{src}V src​ and V tgt V_\textit{tgt}V tgt​ refer to the subset of values that are within the domain of definition of src\textit{src}src and tgt\textit{tgt}tgt respectively. For this chapter, we will assume V src=V tgt=V V_\textit{src} = V_\textit{tgt} = V V src​=V tgt​=V; in other words, minIR graphs are IO-free (this is w.l.o.g., see discussion after Proposition .) and all values are linear1 (this is definitely not w.l.o.g.!).

As a result of this assumption, the subcategory of minIR graphs that we consider forms a rigid category, as introduced by Danos et al. Danos, 2014Vincent Danos, Reiko Heckel and Pawel Sobocinski. 2014. Transformation and Refinement of Rigid Structures. The definition, which we reproduce here, is given in terms of morphisms that intersect all components of the codomain. We refer to Danos, 2014Vincent Danos, Reiko Heckel and Pawel Sobocinski. 2014. Transformation and Refinement of Rigid Structures for the precise definition of that notion: in the context of linear-valued minIR graphs, this is equivalent to requiring that the image of the graph homomorphism intersects every connected component of the codomain.

Definition 4.1Rigid Category

A category C\mathbb CC is rigid if for all morphisms A→h B A \xrightarrow{h} BA h
​B in C\mathbb CC that intersects all components of B B B and for all A→g C A \xrightarrow{g} CA g
​C that factorises as g=f∘h g = f \circ h g=f∘h, then B→f C B \xrightarrow{f} C B f
​C is unique.

In other words, there is a unique way to extend a morphism A→g C A \xrightarrow{g} CA g
​C to a morphism B→f C B \xrightarrow{f} C B f
​C, if it exists. If we interpret A AA and B B B as graph patterns that we are interested in matching with A⊆B A \subseteq BA⊆B, then rigidity guarantees that there is (at most) a unique way to extend a match morphism A→G A \to GA→G into a match morphism on the larger pattern B→G B \to G B→G.

The linearity assumption also has other useful consequences. Every linear value has exactly one use and one definition. As a result, all linear paths are disjoint and form a partition of the values of the graph. They correspond to the paths that form the connected components of the fully split graph, i.e. the graph obtained by splitting every operation. We call the number of linear paths (and hence the number of connected components in the fully split graph) the circuit width, written w i d t h(G)width(G)w idth(G). We also use the linear path decomposition to define circuit depth, written d e p t h(G)depth(G)depth(G), as the longest linear path in G GG.

As discussed in section 3.4, minIR rewrites are instantiated from transformation rules by minIR match morphisms m: P→G m: P \to Gm: P→G. Restricting our considerations to linear-valued minIR graphs has the further implication that all such match morphsisms m mm will be injective. We call m mm an embedding and write it using greek letters and a hooked arrow

m=φ: P↪G. m = \varphi: P \hookrightarrow G.m=φ: P↪G.

Finding such embeddings P↪G P \hookrightarrow G P↪G is the pattern matching problem that we are solving. This problem is equivalent to finding minIR subgraphs H⊆G H \subseteq G H⊆G of G GG such that H H H is isomorphic to the pattern P≃H P \simeq H P≃H.

#### Convexity #

According to Proposition ., a necessary condition for a subgraph H H H to define a valid minIR rewrite is convexity. In this chapter we weaken this requirement and propose a condition based on circuit width:

Proposition 4.1Necessary condition for convexity

Let φ: P↪G\varphi: P \hookrightarrow Gφ: P↪G be an embedding of a pattern P P P into a linear-valued minIR graph G GG such that φ(P)\varphi(P)φ(P) is a convex subgraph of G GG. Then for every subgraph H⊆G H \subseteq G H⊆G such that φ(P)⊆H\varphi(P) \subseteq Hφ(P)⊆H, it holds that w i d t h(P)≤w i d t h(H).width(P) \leq width(H).w idth(P)≤w idth(H).

Up to isomorphism, we can assume P⊆G P \subseteq G P⊆G. Suppose there is H⊆G H \subseteq G H⊆G such that P⊆H P \subseteq H P⊆H and w i d t h(P)>w i d t h(H)width(P) > width(H)w idth(P)>w idth(H). Let L P, L H⊆P(V)\mathcal{L}_P, \mathcal{L}_H \subseteq \mathcal{P}(V)L P​,L H​⊆P(V) be partitions of V P V_P V P​ and V H V_H V H​ into sets of values that are on the same linear path of P P P and H H H respectively. It must hold that for all ℓ∈L P\ell \in \mathcal{L}_Pℓ∈L P​ there is ℓ′∈L H\ell' \in \mathcal{L}_Hℓ′∈L H​ such that ℓ⊆ℓ′\ell \subseteq \ell'ℓ⊆ℓ′, because P⊆H P \subseteq H P⊆H and operation splitting is preserved under embeddings. As the map from L P\mathcal{L}_PL P​ to L H\mathcal{L}_HL H​ cannot be injective, there must be ℓ 1, ℓ 2∈L P\ell_1, \ell_2 \in \mathcal{L}_Pℓ 1​,ℓ 2​∈L P​ and ℓ′∈L H\ell' \in \mathcal{L}_Hℓ′∈L H​, such that ℓ 1≠ℓ 2\ell_1 \neq \ell_2ℓ 1​=ℓ 2​ and ℓ 1, ℓ 2⊆ℓ′\ell_1, \ell_2 \subseteq \ell'ℓ 1​,ℓ 2​⊆ℓ′. We conclude that there must be a path in the fully split graph of H H H between a value of ℓ 1\ell_1ℓ 1​ and a value of ℓ 2\ell_2ℓ 2​ that is not in the fully split graph of P P P. Given that P P P is convex, this path must be in P P P, which contradicts the preservation of operation splitting under embeddings.

In this chapter, whenever we define a subgraph H⊆G H \subseteq G H⊆G of a graph G GG, we will assume that H H H satisfies the above weakened convexity condition.

The converse of Proposition ., however, is not true. The pattern-matching technique presented below will find a strict superset of convex embeddings. To restrict considerations to convex embeddings, it suffices to filter out the non-convex ones in a post-processing step.

#### Ignoring minIR Hierarchy #

So far, we have omitted discussing one part of the minIR structure: the nested hierarchy of operations. Syntactically, the hierarchy formed by p a r e n t parentparent relations between minIR operations can be viewed as just another value type that operations are incident to: parent operations define an additional output that children operations consume as additional input. Because of the bijectivity requirement of minIR morphisms on parent-child relations of Definition ., these parent-child relations behave, in fact, like linear values – and hence do not violate the linearity assumption we have imposed.

However, by treating them as such, we have further weakened the constraints on pattern embeddings. We do not enforce that boundary values must be in the same regions or that parent-child relations cannot be boundary values. Similarly to convexity, we defer checking these properties to a post-processing step.

#### Further assumptions (harmless) #

We will further simplify the problem by making presentation choices that do not imply any loss of generality. First of all, we assume that all patterns have the same width w w w and depth d dd, are connected graphs and have at least 2 operations. These conditions can always be fulfilled by adding “dummy” operations if necessary. Embeddings of disconnected patterns can be computed one connected component at a time.

We will further assume that all operations are on at most two linear paths (and thus in particular, have at most 4 endpoints). Operations on Δ>2\Delta > 2Δ>2 linear paths can always be broken up into a composition of Δ−1\Delta-1Δ−1 operations, each on two linear paths as follows:

Expressing an operation on Δ=3\Delta = 3Δ=3 linear paths as a composition of two operations on 2 linear paths.

This transformation leaves circuit width unchanged but may multiply the graph depth by up to a factor Δ\DeltaΔ.

We furthermore define the set of all port labels

P a l l=⋃ν∈Γ P ν P_{all} = \bigcup_{\nu \in \Gamma} P_\nu P a ll​=ν∈Γ⋃​P ν​

so that we can associate every operation endpoint in a minIR graph G GG with a port label from the set P a l l P_{all}P a ll​. We further endow the labels P a l l P_{all}P a ll​ with a total order (for instance, based on the string index values). The total order on P a l l P_{all}P a ll​ then induces a total order on the paths v 1⋯v k∈V∗v_1\cdots v_k \in V^\ast v 1​⋯v k​∈V∗ in G GG that start in the same value v 1 v_1 v 1​: the paths are equivalently described by the sequence of port labels of the operations traversed. These form strings in P a l l∗P_{all}^\ast P a ll∗​, which we order lexicographically. Given a root value r r r, for every value v v v in G GG there is thus a unique smallest path from r r r to v v v in G GG2. This path is invariant under isomorphism of the underlying graph (i.e. relabelling of the values and operations but preserving the port labels). With this we conclude the discussions of the specificities of minIR graphs related to typing, linearity and hierarchy, and the related assumptions that we are making.

To summarise, minIR graphs as they are considered in this chapter are hypergraphs (Definition 3.1) that satisfy the following properties

- every vertex (value) is incident to exactly two hyperedges (operations). It is the target of one hyperedge (its definition) and the source of another one (its use),
- every hyperedge is incident to at most four vertices,
- every hyperedge can be split in a unique way (and invariant under isomorphism) into at most two split operations, with each at most two endpoints.

When modelling subgraphs of IO-free minIR graphs (typically patterns for pattern matching), some hyperedge connections at the boundary of the subgraph will be missing. We say a value is open if a use or define operation is missing (i.e. it is a boundary value in a minIR subgraph).

We will simplify refer to hypergraphs that satisfy the above assumptions as graphs. In the unique instance of this chapter where a graph that does not satisfy this construction is referred to, we will specifically call it a simple graph.

We conclude with the following notable bound on circuit width.

Proposition 4.2Bound on circuit width

Let G GG be a graph with n odd n_\textrm{odd}n odd​ operations of odd arity (i.e. ∣d e f(o)+u s e(o)∣|def(o) + use(o)|∣de f(o)+use(o)∣ is odd) and n ω n_\omegan ω​ open values. Then, the circuit width of G GG is

w i d t h(G)=⌊(n odd+n ω)/2⌋.width(G) = \lfloor(n_\textrm{odd} + n_\omega) / 2\rfloor. w idth(G)=⌊(n odd​+n ω​)/2⌋.

For any linear path P⊆V∗P \subseteq V^\ast P⊆V∗ in G GG consider its two ends v 1 v_1 v 1​ and v 2 v_2 v 2​, i.e. the two values in P P P with only one neighbouring value in P P P (by definition linear paths cannot be empty). In the fully split graph of G GG, these values are either open or must be connected to two operations. In the latter case, at least one of the operations must have a single endpoint (otherwise by acyclicity, the operation would have two neighbours).

In a fully split graph, operations with a single endpoint result from a split operation with an odd number of endpoints. We conclude that for every linear path, there are either two operations with an odd number of endpoints in G GG, or one such operation and one open value, or two open values. The result follows.

1. This restriction is necessary for our results: copyable values may admit an arbitrary number of adjacent hyperedges. As a result, minIR graph pattern matching with copyable values is a generalisation of the subgraph isomorphism problem, a well-known NP-complete problem Cook, 1971Stephen A. Cook. 1971. The complexity of theorem-proving procedures. In Proceedings of the third annual ACM symposium on Theory of computing - STOC ’71. ACM Press, 151--158. doi: 10.1145/800157.805047. The approach generalises to non-linear types, but the complexity analysis no longer holds (we pay a computational price for every non-linear value matched). ↩︎

2. Remark that the ordering of the operations thus defined is a particular case of a depth-first search (DFS) ordering of the graph: given an operation o oo that has been visited, all its descendants will be visited before proceeding to any other operation. ↩︎

## 4.3. Tree reduction

We reduce the problem of graph pattern matching to matching on rooted trees – as we will see in section 4.6, a much simpler problem to solve. The map between graphs and (rooted) trees is given by rooted dual trees. Call G GG tree-like if G GG is connected and the underlying undirected graph G U G_UG U​ of G GG is acyclic.

Definition 4.2Rooted dual tree

Let G GG be a tree-like graph with operations O O O. Then given a root operation r∈O r \in O r∈O, the rooted dual tree of G GG rooted at r r r, written τ r(G)\tau_r(G)τ r​(G) is the tree given by

- the nodes of the tree are the operations O O O of G GG,
- the parent and children of o∈O o \in Oo∈O are the operations that share a value with o oo in G GG; the parent is the unique operation on the path from o oo to r r r,
- the children of an operation are ordered according to the port labels.

Unlike graphs, tree nodes are identified uniquely by their path from the root. Trees isomorphic as graphs with identical root are thus considered equal.

#### A tree reduction using path splitting #

To reduce a graph G GG to a tree using the rooted dual tree construction, it suffices to reduce G GG to a tree-like graph. The following result shows that this can always be achieved by repeatedly applying operation splitting transformations.

Proposition 4.6Path splitting

A tree-like graph can be obtained from any connected graph G GG by applying operation splittings. The resulting graph is a path-split graph (PSG) of G GG.

Consider the undirected simple graph I\mathcal{I}I, where vertices are linear paths, and there is an edge between two linear paths for every operation that belongs to both paths. We call I\mathcal{I}I the interface graph of G GG.

Splitting an operation o oo in a graph G GG corresponds to removing the corresponding edge in I\mathcal{I}I. On the other hand, the underlying undirected graph G U G_UG U​ of G GG has a cycle if and only if there is a cycle in I\mathcal{I}I. Indeed, a cycle in G U G_UG U​ cannot belong to a single linear path in G GG, by acyclicity of minIR graphs. There is, therefore, a cycle of operations that span multiple linear paths, thus forming a cycle in I\mathcal{I}I.

Hence, the operations to be split to turn G GG into a tree-like graph are given by the set of edges E−E^-E− in I\mathcal{I}I that must be removed to obtain a spanning tree of I.\mathcal{I}.I.1

As we consider typed graph, the splitting of an operation is unique; however the choice of spanning tree of I\mathcal I I is not unique, and thus multiple PSGs exist for a given graph G GG.

If G′G'G′ is a PSG of some graph G GG, then call an operation o oo of G GG an anchor operation if it is on two linear paths and it is not split in G′G'G′. The set of all anchors operations π⊆O\pi \subseteq O π⊆O fully determines the path-split graph. We write G π=G′G^\pi = G'G π=G′ for the PSG of G GG obtained by anchors π\pi π.

Proposition 4.3Rooted dual trees are ternary
Consider a PSG G π G^\piG π of a graph G GG. There is a root operation r∈O r \in O r∈O such that the rooted dual tree τ r(G π)\tau_r(G^\pi)τ r​(G π) is a ternary tree, i.e. every node of τ r(G π)\tau_r(G^\pi)τ r​(G π) has at most three children.

We have assumed in section 4.2 that every operation in G GG is on at most two linear paths and thus can be connected to at most four values. Each value is linear and hence connected to at most one other operation. It results that every operation in τ r(G π)\tau_r(G^\pi)τ r​(G π) has at most four neighbouring operations – one parent and three children. A tree leaf can be chosen as the root operation to ensure the root does not have four children.

We can make the path splitting transformation G→G π G \to G^\piG→G π reversible by separately storing the set of split operations in G π G^\piG π that correspond to a single operation in G GG. As every operation of G GG can get split in at most two split operations, we can store the pairs of split operations in G π G^\piG π that correspond to an operation in G GG in a partial map that defines weights for (a subset of) the operations O π O^\pi O π of G π G^\piG π:

s p l i t: O π⇀P∗.split: O^\pi \rightharpoonup P^\ast.s pl it: O π⇀P∗.

This maps a split operation o oo to the unique undirected path in G π G^\piG π from o oo to the other half of the split operation.

This defines a map σ 1:(G π, s p l i t)↦G\sigma_1: (G^\pi, split) \mapsto G σ 1​:(G π, s pl it)↦G, the inverse of the path splitting transformation G→G π G \to G^\piG→G π.

#### Contracted path-split graphs #

We can further simplify the structure of the data of a PSG by contracting all operations of G π G^\piG π that are on a single linear path. The result is the contracted path-split graph (cPSG) of G GG, written c(G π)c(G^\pi)c(G π).

We employ a similar trick as above to make this transformation reversible, this time by introducing weights on the values of c(G π)c(G^\pi)c(G π) that store the string of operations that were contracted2

c o n t r a c t: V C⇀(Γ π)∗contract: V_C \rightharpoonup (\Gamma^\pi)^\astcont r act: V C​⇀(Γ π)∗

where V C V_C V C​ are the values of c(G π)c(G^\pi)c(G π) and Γ π\Gamma^\piΓ π are the optypes of operations in G π G^\piG π, i.e. the optypes of the minIR graph G GG along with the optypes of the split operations. This defines a second map σ 2:(c(G π),c o n t r a c t)↦G π\sigma_2: (c(G^\pi), contract) \mapsto G^\pi σ 2​:(c(G π),cont r act)↦G π that is the inverse of the path-split graph contraction transformation c(⋅)c(\cdot)c(⋅). In summary, we have the composition

(c(G π),c o n t r a c t, s p l i t)→σ 2×I d(G π, s p l i t)→σ 1 G.(c(G^\pi), contract, split) \xrightarrow{\sigma_2 \times Id} (G^\pi, split) \xrightarrow{\sigma_1} G.(c(G π),cont r act, s pl it)σ 2​×I d
​(G π, s pl it)σ 1​
​G.

Contracted PSGs are particularly useful for the study of the asymptotic complexity of the pattern matching algorithm we propose, as they have a very regular structure. This is expressed by the following proposition that further extends the statement of Proposition 4.3:

Proposition 4.4Contracted PSG
Consider a PSG G π G^\piG π of a graph G GG. There is a root operation r∈O r \in O r∈O such that the rooted dual tree of the contracted PSG τ r(c(G π))\tau_r(c(G^\pi))τ r​(c(G π)) is a ternary tree with w i d t h(G)−1 width(G) - 1 w idth(G)−1 nodes.

That the tree is ternary follows from Proposition 4.3. Every node of the tree corresponds to an operation in c(G π)c(G^\pi)c(G π), which is on exactly two linear paths. As a result of acyclicity of the tree, a tree of k k k nodes spans k+1 k+1 k+1 linear paths – and hence, we conclude k=w i d t h(G)−1 k = width(G) - 1 k=w idth(G)−1.

We conclude the construction presented in this section with the following result, expressing graph pattern matching in terms of tree equality:

Proposition 4.5Reduction to Tree Pattern matching

Let P P P be a pattern graph and G GG a graph. Let P π P^\pi P π be a PSG of G GG. There is an embedding P↪G P \hookrightarrow G P↪G if and only if there is H⊆G H \subseteq G H⊆G and a PSG H π′H^{\pi'}H π′ of H H H such that

τ(c(P π))=τ(c(H π′))\tau(c(P^\pi)) = \tau(c(H^{\pi'}))τ(c(P π))=τ(c(H π′))

and the trees have equal weight maps s p l i t splits pl it and c o n t r a c t contractcont r act.

The proof of this follows directly from our construction, the unicity of trees under isomorphism and the bijection between the graphs P, H P, H P, H and their cPSGs.

We have thus successfully reduced the problem of pattern matching to the problem of matching on trees. Given that the ordering of children of a node in a tree is fixed, checking trees for equality is a simple matter of checking node and weight equality, one node (and edge) at a time.

We conclude this section with a figure summarising the constructions we have presented.

A graph G GG, along with the path-split graph G π G^\piG π, the contracted path-split graph c(G π)c(G^\pi)c(G π) and their rooted dual trees. The anchor operations are d dd (grey) and e ee (red). The root of the rooted dual trees is e ee.

1. It is a simple result from graph theory that such a set of edges always exists – it suffices to remove one edge from every cycle in the graph. ↩︎

2. Because all contracted operations apply on a single, shared, linear path, they indeed form a string of operations. ↩︎

## 4.4. Canonicalising the tree reduction

The reduction of graph matching to ternary trees from the previous section is a big step towards an algorithm for graph matching. However, Proposition 4.5 is expressed in terms of existence of PSGs – it is as yet unclear how the trees can be constructed. This is the purpose of this section.

We introduce for this purpose a canonical, that is, invariant under isomorphism, choice of PSG G π G^\piG π of G GG. The result is a unique canonical transformation G↦G π↦c(G π)G \mapsto G^\pi \mapsto c(G^\pi)G↦G π↦c(G π) from G GG to a cPSG that we can use for pattern matching.

We proceed by using the total order that we have defined on port labels and can be extended lexicographically to paths outgoing from a shared root operation (see section 4.2 for more details). Whenever more than one path from r r r to o oo exist in G GG, it suffices to consider the smallest one. For a choice of root operation r∈O, r \in O, r∈O, we thus obtain a total order of all operations O O O in G GG.

We then restrict our attention to operations on two linear paths and consider them in order. We keep track of linear paths that have been visited and proceed as follows to determine whether o∈O o \in Oo∈O must be split:

- if o oo is on a linear path that was not seen before, it is left unchanged and the set of visited linear paths is updated;
- otherwise, i.e. o oo is on two linear paths that have already been visited, the operation is split, resulting in two operations on a single linear path.

The pseudocode `CanonicalPathSplit` implements this algorithm. We use `Operations(G)` to retrieve all the operations on the graph `G` and `LinearPaths(G, op)` to retrieve the linear paths of the operation `op`. The linear paths are identified using integer indices that can be pre-computed and stored in linear time in the graph size. `SplitOperation(G, op)` returns the graph resulting from splitting `op` into two operations on a single linear path. Finally, `PathAsPortLabels(G, root, v)` returns the string of the port labels that encode the path from `root` to `v` in the graph `G`. The strings are ordered lexicographically. The non-capitalized functions `set`, `union`, `sort`1, `len`, and `issubset` have their standard meanings.

     1def CanonicalPathSplit(G: Graph, root: Operation) -> Graph:
     2  new_G := G
     3  all_operations := Operations(G)
     4  sorted_operations := sort(
     5      all_operations,
     6      sort_key= lambda v: PathAsPortLabels(G, root, v)
     7  )
     8
     9  # keep track of the visited linear paths
    10  seen_paths := set()
    11  for op in sorted_operations:
    12    # Get the (pre-computed) indices of the linear paths
    13    op_linear_paths := LinearPaths(G, op)
    14    if len(op_linear_paths) == 2:
    15      if issubset(op_linear_paths, seen_paths):
    16        # The two linear paths of `op` are already visited
    17        new_G = SplitOperation(new_G, op)
    18      else:
    19        # Mark the new linear paths as visited
    20        seen_paths = union(seen_paths, op_linear_paths)
    21  return new_G

The following figure shows an example of splitting a graph into its canonical PSG using `CanonicalPathSplit`.

Splitting a graph into its canonical PSG. Ports are ordered counter-clockwise on each edge, and numbered according to the lexicographic order of the paths from root to the port, as returned by `PathAsPortLabels`. This induces an order on the hyperedges, reflected in the alphabetic order of the edge labels. Linear paths are formed by ports in a horizontal line (as marked by the dotted lines). Vertex root is chosen as the root of the canoncal splitting. Vertices d and g are not split because they are the smallest edges that contain the fourth, respectively first linear path.

Proposition 4.7Correctness of `CanonicalPathSplit`
For a graph G GG, the graph returned by `CanonicalPathSplit(G)` is a valid PSG of G. G.G. It is deterministic and invariant under isomorphism of G GG. The runtime of `CanonicalPathSplit` is O(∣G∣)O(|G|)O(∣G∣), where ∣G∣|G|∣G∣ is the number of operations in the graph G GG.

Let G π G^\piG π be the graph returned by `CanonicalPathSplit(G)`. From the discussion in the proof of Proposition 4.6, we know it is sufficient to show that the interaction graph I\mathcal{I}I of G π G^\piG π is acyclic and connected.

I\mathcal{I}I is acyclic. If there was a cycle in I\mathcal{I}I, then there would be operations o 0,…,o k−1 o_0, \dots, o_{k-1}o 0​,…,o k−1​ in G GG that pairwise (o i, o i+1 m o d k)(o_i, o_{i+1\, mod\, k})(o i​,o i+1 mod k​) share a linear path. One of these operations must be considered last in the `for` loop of lines 11–20, suppose it is o k−1 o_{k-1}o k−1​. But every linear path of o k−1 o_{k-1}o k−1​ is either also a linear path of o k−2 o_{k-2}o k−2​ or a linear path of o 1 o_{1}o 1​: o k−1 o_{k-1}o k−1​ thus does not satisfy the condition on line 15, and thus cannot be in I\mathcal{I}I, a contradiction. Hence I\mathcal{I}I is acyclic.

I\mathcal{I}I is connected. We proceed inductively to show the following invariant for the main `for` loop (lines 11–20): for all linear paths in `seen_paths`, there is a path in I\mathcal{I}I to a linear path of the root operation. `seen_paths` is only modified on line 20. If `op` is the root operation, then trivially there is a path from the linear paths `op_linear_paths` to a linear path of the root operation. Otherwise, we claim that there must be one of the paths in `op_linear_paths` that is already in `seen_paths`. From there it follows that there is a path in I\mathcal{I}I from the root path to the unseen linear path, given by the path to the linear path in `seen_path` followed by the edge in I\mathcal{I}I that corresponds to `op`.

By connectedness of G GG, there is a path from the root operation to `op`. The path is not empty because `op` is not the root operation, so we can consider the prefix of the path of all operations excluding `op`. Call `op'` the last operation preceding `op` and `op_linear_paths'` its linear paths. Two successive operations on a path must share a linear path: `op_linear_paths` ∩\cap∩ `op_linear_paths'` cannot be empty. According to line 4, `op'` must have been visited before `op`, thus `op_linear_paths'` ⊆\subseteq⊆ `seen_paths`. It follows that at least one element of `op_linear_paths` must be in `seen_paths`.

Determinstic and isomorphism invariant. The pseudocode above is deterministic and only depends on paths in G GG encoded as strings of port labels, which are invariant under isomorphism.

Runtime complexity. Lines 2 and 3 run in O(∣G∣)O(|G|)O(∣G∣) time. With the exception of the `sort` function on lines 4–7, every other line can be run in O(1)O(1)O(1) time:

- lines 13 and 15 run in constant time because the size of `op_linear_paths` is always at most 2;
- line 20 (and the `in` check on line 15) can be run in constant time by representing the `seen_paths` set as a fixed-size boolean array of size w w w, with the i ii-th bit indicating whether the i ii-th linear path has been seen;
- line 17 is a constant time transformation if we allow in-place modification of `new_G`.

The `for` loop will run ∣G∣|G|∣G∣ iterations, for a total of O(∣G∣)O(|G|)O(∣G∣) runtime. Finally, the sorting operation would naively take time O(∣G∣log⁡∣G∣)O(|G| \log |G|)O(∣G∣log∣G∣). However, given that the ordering is obtained lexicographically from the paths starting at the root, we can obtain the sorted list of operations by depth-first traversal of the graph starting at the root. The result follows.

Using `CanonicalPathSplit`, we can now sketch what the pattern matching algorithm should look like. For each pattern, we first compute their canonical PSG for an arbirary choice of pattern root operation; then, given a graph G GG, we can find all embeddings of patterns into G GG by iterating over all possible PSGs within G GG. Naively, this involves enumerating all posible subgraphs of G GG, and then for each of them, iterating over all possible root choices.

This can be significantly sped up by realising that many of the PSGs that are computed when iterating over all possible subgraphs and root choices are redundant2. We will see in the next section that we can i) iterate once over all possible root choices r r r in G GG and ii) introduce a new procedure `AllPathSplits` that will efficiently enumerate all possible rooted ual trees of PSGs that are rooted in r r r for subgraphs within G GG. In the process, we will also see that we can replace the tree equality check of line 12 with a subtree inclusion check, further reducing the number of PSGs that must be considered.

Naive pattern matching.

     1# Precompute all PSGs
     2allT = [CanonicalPathSplit(
     3    P, root_P
     4) for (P, root_P) in patterns]
     5
     6for S in Subgraphs(G):
     7  for root_S in Operations(S):
     8    TG = CanonicalPathSplit(
     9        S, root_S
    10    )
    11    for T in allT:
    12      if T == TG:
    13        yield T

Improved using `AllPathSplits` (section 4.5).

     1# Precompute all PSGs
     2allT = [CanonicalPathSplit(
     3    P, root_P
     4) for (P, root_P) in patterns]
     5
     6for root_G in Operations(G):
     7  for TG in AllPathSplits(
     8      G, root_G
     9  ):
    10    for T in allT
    11      # Replace == with subtree
    12      if IsSubTree(T, TG)
    13        yield T

1. The `sort_key` parameter of the `sort` function defines the total order according to which the elements are sorted, from smallest to largest. ↩︎

2. Think for example of the same root operation r r r that is considered repeatedly for every overlapping subgraph of G GG that contains r r r. ↩︎

## 4.5. Enumerating all path-split graphs

The `CanonicalPathSplit` procedure in the previous section defines for all graphs G GG and choice of root operation r r r a canonical PSG G π G^\piG π, and thus a canonical set of anchors π\pi π that we write as π r(G)=π.\pi_r(G) = \pi. π r​(G)=π.

Instead of `CanonicalPathSplit`, we can equivalently consider a `CanonicalAnchors` procedure, which computes π r(G)\pi_r(G)π r​(G) directly instead of the graph G π r(G)G^{\pi_r(G)}G π r​(G).

We formulate this computation below, using recursion instead of a `for` loop. This form generalises better to the `AllAnchors` procedure that we will introduce next.

The equivalence of the `CanonicalAnchors` procedure with `CanonicalPathSplit` follows from the observation made in section 4.2 that ordering operations in lexicographic order of the port labels is equivalent to a depth-first traversal of the graph.

`CanonicalAnchors` implements a recursive depth-first traversal (DFS), with the twist that the recursion is explicit only on the anchor nodes and otherwise relying on the lexicographic ordering just like in `CanonicalPathSplit`: lines 5–15 of `CanonicalAnchors` correspond to the iterations of the `for` loop (line 11–20) of `CanonicalPathSplit` until an anchor operation is found (i.e. the `else` branch on lines 18–20 is executed). From there, the graph traversal proceeds recursively.

We introduce the `ConnectedComponent`, `Neighbours` and `RemoveOperation` procedures; the first returns the connected component of the current operation, whereas the other two procedures are used to traverse, respectively modify, the graph G GG. Importantly, `Neighbours(G, op)` returns the neighbours of `op` ordered by port label order.

To ensure that the recursive DFS does not visit the same operation twice, we modify the graph with `RemoveOperation` on lines 11 and 15, ensuring that no visited operation remains in `G`. As a consequence, `CanonicalAnchors` may be called on disconnected graphs, which explains why an additional call to `ConnectedComponent` (line 4) is required.

Proposition 4.8Equivalence of `CanonicalPathSplit` and `CanonicalAnchors`

Let G GG be a connected graph and let r r r be a root operation in G GG. Then `CanonicalAnchors` maps the graph to the canonical anchor set:

(G, r,{})↦(π(G)r, L,∅),(G, r, \{\}) \mapsto (\pi(G)_r, L, \varnothing),(G, r,{})↦(π(G)r​,L,∅),

where L LL is the set of all paths in G GG and ∅\varnothing∅ designates the empty graph.

The proof follows directly from the previous paragraphs.

     1def CanonicalAnchors(
     2    G: Graph, root: Operation, seen_paths: Set[int]
     3) -> (Set[Operation], Set[int], Graph):
     4  operations = Operations(ConnectedComponent(G, root))
     5  # sort by PathAsPortLabels, as previously
     6  sorted_operations := sort(operations)
     7  operations_queue := queue(sorted_operations)
     8
     9  # Skip all operations that are not anchors
    10  op := operations_queue.pop() # never emtpy, contains root
    11  G = RemoveOperation(G, op)
    12  while len(LinearPaths(G, op)) == 1 or
    13        issubset(LinearPaths(G, op), seen_paths):
    14    op = operations_queue.pop() or return ({}, {}, G)
    15    G = RemoveOperation(G, op)
    16
    17  # op is anchor, update seen_paths and recurse
    18  seen_paths = union(seen_paths, LinearPaths(G, op))
    19  anchors := [op]
    20  # sort by port labels
    21  for child in Neighbours(G, op):
    22    (new_anchors, seen_paths, G) = CanonicalAnchors(
    23        G, child, seen_paths
    24    )
    25    anchors += new_anchors
    26
    27  return (anchors, seen_paths, G)

#### Maximal PSGs #

In addition to “simplifying” the data required to define path splitting, the definition of PSGs using anchor operations has another advantage that is fundamental to the pattern matching algorithm.

Consider the rooted dual tree τ r(G π)\tau_r(G^\pi)τ r​(G π) of a PSG with root operation r r r in G π G^\piG π. Recall that tree nodes are uniquely identified by their path from the root and thus are considered equal if they are isomorphic as graphs. We can in the same way define a tree inclusion relation ⊆\subseteq⊆ on rooted dual trees that corresponds to checking that the trees have the same root and that the left-hand side is isomorphic to a subtree of the right-hand side. We also require that the operation weights given by the s p l i t splits pl it map s p l i t splits pl it map coincide on the common subtree.

Proposition 4.9Maximal PSG

Let G GG be a connected graph, π\pi π a set of operations in G GG and r∈π r \in \pi r∈π a root operation. Consider the set G π={H⊆G∣π r(H)=π}.\mathcal{G}_\pi = \{H \subseteq G \mid \pi_r(H) = \pi \}.G π​={H⊆G∣π r​(H)=π}.

There is a subgraph M⊆G M \subseteq G M⊆G such that for all subgraphs H∈G X H \in \mathcal{G}_X H∈G X​: H⊆M H \subseteq M H⊆M. Furthermore, for all graph P P P, there is r′r'r′ and π′=π(P)r′\pi' = \pi(P)_{r'}π′=π(P)r′​ such that

P≃H∈G X⇔τ r′(P π′)⊆τ r(M π).P \simeq H \in \mathcal{G}_X \quad\Leftrightarrow\quad \tau_{r'}(P^{\pi'}) \subseteq \tau_r(M^\pi).P≃H∈G X​⇔τ r′​(P π′)⊆τ r​(M π).

We call M π M^\pi M π the maximal PSG with anchors π\pi π in G GG.

The proof gives an explicit construction for M M M.

Assume G X≠∅\mathcal{G}_X \neq \varnothing G X​=∅, otherwise the statement is trivial.

Construction of M M M. Let L LL be the set of linear paths in G GG that go through at least one operation in π\pi π. Consider the set of operations O L O_L O L​ in G GG given by the operations whose linear paths are contained in L LL. This defines a subgraph G∣O L G|_{O_L}G∣O L​​ of G GG. Since G X≠∅\mathcal{G}_X \neq \varnothing G X​=∅, there exists H∈G X H \in \mathcal{G}_X H∈G X​. By assumption, H H H is connected, and thus the anchors π\pi π of H H H are connected in H H H. There is therefore a connected component M⊆G∣O L M \subseteq G|_{O_L}M⊆G∣O L​​ that contains the set π\pi π.

Well-definedness of M π M^\pi M π. Consider the PSG M π M^\pi M π of M M M. We must show that M π M^\pi M π is a tree-like graph for the proposition statement to be well-defined. In other words, we must show that the interaction graph I\mathcal{I}I of M π M^\pi M π is acyclic and connected. M M M is connected by construction, which implies connectedness of M π M^\pi M π and thus of I\mathcal{I}I. It is acyclic because w i d t h(M)=∣π∣+1 width(M) = |\pi| + 1 w idth(M)=∣π∣+1 and M M M has exactly ∣π∣|\pi|∣π∣ operations on more than one linear path. I\mathcal{I}I is a thus a tree.

H⊆M H \subseteq M H⊆M. For any subgraph H∈G X H \in \mathcal{G}_X H∈G X​, its operations must be contained in O L O_L O L​. Since any H∈G H \in \mathcal{G}H∈G is connected and contains π\pi π, it must further hold that H⊆M H \subseteq M H⊆M.

We can now prove the ⇔\Leftrightarrow⇔ equivalence of (1).

⇐\Leftarrow⇐: If τ r′(P π′)⊆τ r(M π)\tau_{r'}(P^{\pi'}) \subseteq \tau_r(M^\pi)τ r′​(P π′)⊆τ r​(M π), then there exists H′⊆M π H' \subseteq M^\pi H′⊆M π with rooted dual tree

τ r(H π)=τ r(P π′).\tau_r(H^\pi) = \tau_r(P^{\pi'}).τ r​(H π)=τ r​(P π′).

Furthermore, by definition of ⊆\subseteq⊆ on rooted trees, a s p l i t splits pl it map is defined on H′H'H′, given by the s p l i t splits pl it map of M π M^\pi M π on the domain H′H'H′. Recall from section 4.3 that there is a map σ\sigma σ that maps

(H′,s p l i t)⟼σ H and(M π, s p l i t)⟼σ M.(H', split) \overset{\sigma}{\longmapsto} H\quad\textrm{and}\quad(M^\pi, split) \overset{\sigma}{\longmapsto} M.(H′,s pl it)⟼σ​H and(M π, s pl it)⟼σ​M.
It merges split operations pairwise, and thus it is immediate that H′⊆M π H' \subseteq M^\pi H′⊆M π implies H⊆M H \subseteq M H⊆M. Thus H∈G X H \in \mathcal{G}_X H∈G X​ and H′=H π H' = H^\pi H′=H π. By construction, one can also derive that P≃H P \simeq H P≃H. The statement follows.

⇒\Rightarrow⇒: Since H∈G X H \in \mathcal{G}_X H∈G X​, we know from point 1 that H⊆M H \subseteq M H⊆M. Thus we can define an injective embedding φ: P→M\varphi: P \to Mφ: P→M.

Operation splitting leaves the set of values from H H H to H π H^\pi H π, as well as from M M M to M π M^\pi M π unchanged. Similarly, there is a bijection between values in H π H^\pi H π and M π M^\pi M π and thus between edges in τ r(H π)\tau_r(H^\pi)τ r​(H π) and τ r(M π)\tau_r(M^\pi)τ r​(M π). The pattern embedding φ\varphiφ hence defines an injective map ϕ E\phi_Eϕ E​ from tree edges in τ r(H π)\tau_r(H^\pi)τ r​(H π) to tree edges in τ r(M π)\tau_r(M^\pi)τ r​(M π). We extend this map to a map on the trees ϕ: τ r(H π)→τ r(M π)\phi: \tau_r(H^\pi) \to \tau_r(M^\pi)ϕ: τ r​(H π)→τ r​(M π) by induction over the nodes set of τ r(H π)\tau_r(H^\pi)τ r​(H π). We start by the root map ϕ(r)=r\phi(r) = rϕ(r)=r. Using ϕ E\phi_Eϕ E​, we can then uniquely define the image of any child node of r r r in τ r(H π)\tau_r(H^\pi)τ r​(H π), and so forth inductively.

We show now that the map ϕ\phiϕ thus defined is injective. Suppose v, v′v, v'v, v′ are nodes in τ r(H π)\tau_r(H^\pi)τ r​(H π) such that ϕ(v)=ϕ(v′)\phi(v) = \phi(v')ϕ(v)=ϕ(v′). By the inductive construction there are paths from the root r r r to v v v and v′v'v′ respectively such that their image under ϕ E\phi_Eϕ E​ are two paths from r r r to ϕ(v)=ϕ(v′)\phi(v) = \phi(v')ϕ(v)=ϕ(v′). But τ r(M π)\tau_r(M^\pi)τ r​(M π) is a tree, so both paths must be equal. By bijectivity of ϕ E\phi_Eϕ E​, it follows v=v′v = v'v=v′, and thus ϕ\phiϕ is injective. Finally, the value and operation weights are invariant under pattern embedding and thus are preserved by definition.

This result means that instead of listing all PSGs for every possible subgraph of G GG, it is sufficient to proceed as follows:

1. for every pattern P P P, fix a root operation r P r_P r P​ and construct the rooted tree dual of the canonical PSG
   τ P:=τ r P(P π r P(P)).\tau_P := \tau_{r_P}(P^{\pi_{r_P}(P)}).τ P​:=τ r P​​(P π r P​​(P)).
2. enumerate every possible root operation r r r in G GG,
3. enumerate every possible sets of anchors π\pi π in G GG with root r r r,
4. for each set π\pi π, find the maximal PSG M M M with anchors π\pi π in G GG, and take its rooted tree dual τ M:=τ r(M π)\tau_M := \tau_r(M^\pi)τ M​:=τ r​(M π),
5. find all patterns P P P such that τ P⊆τ M\tau_P \subseteq \tau_M τ P​⊆τ M​.

In other words, if `AllAnchors` is a procedure that enumerates all possible sets of anchors π\pi π in G GG and `MaximalPathSplit` computes the maximal PSG M M M as presented in the proof of Proposition 4.9, then `AllPathSplits(G)` can simply be obtained by calling `AllAnchors` and then returning their respective maximal PSGs in G GG:

    def AllPathSplits(G: Graph, root: Operation) -> Set[Graph]:
      all_anchors = AllAnchors(G, root)
      return {MaximalPathSplit(G, pi) for pi in all_anchors}

#### The missing piece: `AllAnchors` #

We can now complete the definition of `AllPathSplits` by defining the `AllAnchors` procedure, which enumerates all possible sets of anchors in G GG given a root operation r r r.

The procedure is similar to `CanonicalAnchors`, described in detail in the previous paragraphs. In addition to the arguments of `CanonicalAnchors`, `AllAnchors` requires a width w≥1 w \geq 1 w≥1 argument. It then returns all sets of at most w w w operations1 that form the canonical anchors of some width-w w w subgraph of G GG with root r r r. The main difference between `CanonicalAnchors` and `AllAnchors` is that the successive recursive calls (line 22 in `CanonicalAnchors`) are replaced by a series of nested loops (lines 42–48 in `AllAnchors`) that exhaustively iterate over the possible outcomes for different subgraphs of G GG. The results of every possible combination of recursive calls are then collected into a list of anchor sets, which is returned.

The part of the pseudocode that is without comments is unchanged from `CanonicalAnchors`. Using Proposition 4.3, we know that we can assume that every operation has at most 3 children, and thus 3 neighbours in `G`, given that the operations equivalent to parent nodes were removed.

     1def AllAnchors(
     2    G: Graph, root: Operation, w: int,
     3    seen_paths: Set[int] = {}
     4) -> List[(Set[Operation], Set[int], Graph)]:
     5  # Base case: return one empty anchor list
     6  if w == 0:
     7    return [({}, {}, G)]
     8
     9  operations = Operations(ConnectedComponent(G, root))
    10  sorted_operations := sort(operations)
    11  operations_queue := queue(sorted_operations)
    12
    13  op := operations_queue.pop()
    14  G = RemoveOperation(G, op)
    15  while len(LinearPaths(G, op)) == 1 or
    16        issubset(LinearPaths(G, op), seen_paths):
    17    op = operations_queue.pop() or return [({}, {}, G)]
    18    G = RemoveOperation(G, op)
    19
    20  seen0 = union(seen_paths, LinearPaths(G, op))
    21  # There are always at most three neighbours: we
    22  # unroll the for loop of CanonicalAnchors.
    23  [child1, child2, child3] = Neighbours(G, op)
    24  # Iterate over all ways to split w-1 anchors over
    25  # the three children and solve recursively
    26  all_anchors = []
    27  for 0 <= w1, w2, w3 < w with w1 + w2 + w3 == w - 1:
    28    for (anchors1, seen1, G1) in
    29        AllAnchors(G, child1, w1, seen0):
    30      for (anchors2, seen2, G2) in
    31          AllAnchors(G1, child2, w2, seen1):
    32        for (anchors3, seen3, G3) in
    33            AllAnchors(G2, child3, w3, seen2):
    34          # Concatenate new anchor with anchors from all paths
    35          anchors = union([op], anchors1, anchors2, anchors3)
    36          all_anchors.push((anchors, seen3, G3))
    37  return all_anchors

We can represent the sequence of recursive calls to `AllAnchors` as a tree. The call tree for the graph used as example to illustrate `CanonicalAnchors` earlier is given on the next page.

We now show correctness of the procedure. Let us write Π r w(G)\Pi_r^w(G)Π r w​(G) for the set of sets of anchors returned by `AllAnchors(G, r, w, {})`.

Proposition 4.11Correctness of AllAnchors

Let G GG be a graph and H⊆G H \subseteq G H⊆G be a subgraph of G GG of width w w w. Let r r r be a choice of root operation in H H H. We have π r(H)⊆Π r w(G).\pi_r(H) \subseteq \Pi_r^w(G).π r​(H)⊆Π r w​(G).

A call tree for an execution of `AllAnchors` on the example graph of the previous figure with w=3 w = 3 w=3. Starting from the root, each node in the tree corresponds to either picking an operation as anchor or not (thus splitting it). Edges are labelled by the values assigned to w w w for the respective children of the source node. One path from root to leaf leads to no solution (it is impossible to find an unseen linear path from operation b bb. The other paths each lead to a valid set of three anchors.

The proof is by induction over the width w w w of the subgraph H H H. The idea is to map every recursive call in `CanonicalAnchors` to one of the calls to `AllAnchors` on lines 29, 31 or 33. All recursive results are concatenated on line 36, and thus, the value returned by `CanonicalAnchors` will be one of the anchor sets in the list returned by `AllAnchors`.

Let H⊆G H \subseteq G H⊆G be a connected subgraph of G GG of width w w w. We prove inductively over w w w that if (X, S′,H′)=(X, S', H') = (X, S′,H′)=`CanonicalAnchors`$(H, r,S)$ then there is a graph G′G'G′ such that H′⊆G′⊆G H' \subseteq G' \subseteq G H′⊆G′⊆G such that

(X, S′,G′)∈(X, S', G') \in(X, S′,G′)∈ `AllAnchors`(G, r, w, S)(G, r, w, S)(G, r, w, S)

for all valid root operations r r r of H H H and all subsets of the linear paths of H H H in `seen_paths`. The statement in the proposition directly follows this claim.

For the base case w=1 w = 1 w=1, `CanonicalAnchors` will return the anchors `anchors = [op]` as defined on line 19: there is only one linear path, and it is already in `seen_paths`, thus for every recursive call to `CanonicalAnchors`, the `while` condition on line 12 will always be satisfied until all operations have been exhausted and empty sets are returned. In `AllAnchors`, on the other hand, The only values of `w1`, `w2` and `w3` that satisfy the loop condition on line 27 for w=1 w = 1 w=1 are `w1` === `w2` === `w3` =0= 0=0. As a result, given the `w` =0=0=0 base case on lines 6–7, the lines 35 and 36 of `AllAnchors` are only executed once, and the definition of `anchors` on line 36 is equivalent to its definition in `CanonicalAnchors`.

We now prove the claim for w>1 w > 1 w>1 by induction. As documented in `AllAnchors`, we can assume that every operation has at most 3 children. This simplifies the loop on lines 21–25 of `CanonicalAnchors` to, at most, three calls to `CanonicalAnchors`.

Consider a call to `CanonicalAnchors` for a graph H⊆G H \subseteq G H⊆G, a root operation r r r in H H H and a set S S S of linear paths. Let w a w_a w a​, w b w_b w b​ and w c w_c w c​ be the length of the values returned by the three recursive calls to `CanonicalAnchors` of line 22 for the execution of `CanonicalAnchors` with arguments H H H, r r r and S S S. Let c a, c b c_a, c_bc a​,c b​ and c c c_cc c​ be the three neighbours of r r r in H H H. If the child c x c_xc x​ does not exist, then one can set w x=0 w_x = 0 w x​=0 and it can be ignored – the argument below still holds in that case. The definition of `seen0` on line 20 in `AllAnchors` coincides with the update to the variable `seen_paths` on line 18 of `CanonicalAnchors`; similarly, the updates to `G` on lines 14 and 18 of `AllAnchors` are identical to the lines 11 and 15 of `CanonicalAnchors` that update `H`. Let the updated `seen_paths` be the set S a S_a S a​, the updated `G` be G a G_aG a​ and the updated H H H be H a H_a H a​, with H a⊆G a H_a \subseteq G_a H a​⊆G a​.

As every anchor operation reduces the number of unseen linear paths by exactly one (using the simplifying assumptions of section 4.2), it must hold that w a+w b+w c+1=w w_a + w_b + w_c + 1 = w w a​+w b​+w c​+1=w. Thus, for a call to `AllAnchors` with the arguments G GG, r r r, w w w and S S S, there is an iteration of the `for` loop on line 27 of `AllAnchors` such that `w1` =w a= w_a=w a​, `w2` =w b= w_b=w b​ and `w3` =w c= w_c=w c​. It follows that on line 29 of `AllAnchors`, the procedure is called recursively with arguments (G a, c a, w a, S a)(G_a, c_a, w_a, S_a)(G a​,c a​,w a​,S a​). From the induction hypothesis, we obtain that there is an iteration of the `for` loop on line 29 in which the values of `anchors1` and `seen1` coincide with the values of the `new_anchors` and `seen_paths` variables after the first iteration of the `for` loop on line 21 of `CanonicalAnchors`. Call the value of `seen1` (and `seen_paths`) S b S_b S b​. Similarly, call the updated value of `G` in `AllAnchors` G b G_bG b​ and the updated value of `G` in `CanonicalAnchors` H b H_b H b​. We have, by the induction hypothesis, that H b⊆G b H_b \subseteq G_b H b​⊆G b​.

Repeating the argument, we obtain that there are iterations of the `for` loops on lines 30 and 32 of `AllAnchors` that correspond to the second and third recursive calls to `CanonicalAnchors` on line 22 of the procedure. Finally, the concatenation of anchor lists on line 36 of `AllAnchors` is equivalent to the repeated concatenations on line 25 of `CanonicalAnchors`, and so we conclude that the induction hypothesis holds for w w w.

We will see that the overall runtime complexity of `AllAnchors` can be easily derived from a bound on the size of the returned list. For this, we use the following result:

Proposition 4.10Number of anchor sets in `AllAnchors`
For a graph G GG, a root operation r r r in G GG and 1≤w≤w i d t h(G)1 \leq w \leq width(G)1≤w≤w idth(G), the length of the list `AllAnchors`(G, r, w)(G, r, w)(G, r, w) is in O(c w⋅w−3/2)O(c^w \cdot w^{-3/2})O(c w⋅w−3/2), where c=6.75 c = 6.75c=6.75 is a constant.

Let C w C_w C w​ be an upper bound for the length of the list returned by a call to `AllAnchors` for width w w w. For the base case w=0 w = 0 w=0, C 0=1 C_0 = 1 C 0​=1. The returned `all_anchors` list is obtained by pushing anchor lists one by one on line 36. We can count the number of times this line is executed by multiplying the length of the lists returned by the recursive calls on lines 28–32, giving us the recursion relation

C w≤∑0≤w 1, w 2, w 3<w w 1+w 2+w 3=w−1 C w 1⋅C w 2⋅C w 3. C_w \leq \sum_{\substack{0 \leq w_1, w_2, w_3 < w\\w_1 + w_2 + w_3 = w - 1}} C_{w_1} \cdot C_{w_2} \cdot C_{w_3}.C w​≤0≤w 1​,w 2​,w 3​<w w 1​+w 2​+w 3​=w−1​∑​C w 1​​⋅C w 2​​⋅C w 3​​.

Since C w C_w C w​ is meant to be an upper bound, we replace ≤\leq≤ with equality above to obtain a recurrence relation for C w C_w C w​. This recurrence relation is a generalisation of the well-known Catalan numbers Stanley, 2015Richard P. Stanley. 2015. Catalan Numbers. Cambridge University Press. doi: 10.1017/CBO9781139871495, equivalent to counting the number of ternary trees with w w w internal nodes: a ternary tree with w≥1 w \geq 1 w≥1 internal nodes is made of a root along with three subtrees with w 1, w 2 w_1,w_2 w 1​,w 2​ and w 3 w_3 w 3​ internal nodes respectively, with w 1+w 2+w 3=w−1 w_1 + w_2 + w_3 = w-1 w 1​+w 2​+w 3​=w−1. A closed form solution to this problem can be found in Aval, 2008Jean-Christophe Aval. 2008. Multivariate Fuss–Catalan numbers. Discrete Mathematics 308, 20 (October 2008, 4660–4669). doi: 10.1016/j.disc.2007.08.100​:

C w=(3 w w)2 w+1=Θ(c w w 3/2)C_w = \frac{{3w \choose w}}{2w + 1} = \Theta \left(\frac{c^w}{w^{3/2}} \right)C w​=2 w+1(w 3 w​)​=Θ(w 3/2 c w​)

satisfying the above recurrence relation with equality, where c=27/4=6.75 c = 27/4 = 6.75c=27/4=6.75 is a constant obtained from the Stirling approximation:

(3 w w)=(3 w)!(2 w)!w!=Θ(1 w)((3 w)3 e 3)w(e 2(2 w)2)w(e w)w=Θ((27/4)w w 1/2).\begin{aligned}{3w \choose w} = \frac{(3w)!}{(2w)!w!} &= \Theta\left(\frac{1}{\sqrt{w}}\right) \Big(\frac{(3w)^3}{e^3}\Big)^{w}\Big(\frac{e^2}{(2w)^2}\Big)^{w}\Big(\frac{e}{w}\Big)^{w}\\ &= \Theta\left(\frac{(27/4)^w}{w^{1/2}}\right).\end{aligned}(w 3 w​)=(2 w)!w!(3 w)!​​=Θ(w
​1​)(e 3(3 w)3​)w((2 w)2 e 2​)w(w e​)w=Θ(w 1/2(27/4)w​).​

To obtain a runtime bound for `AllAnchors`, it is useful to identify how much of G GG needs to be traversed. If we suppose all patterns have at most depth d dd, then it immediately follows that any operation in G GG that is in the image of a pattern embedding must be at most a distance d dd away from an anchor operation. We can thus equivalently call `AllAnchors` on a subgraph of G GG such that no linear path is longer than 2 d 2d2d. We thus obtain the following runtime.

Proposition 4.12Runtime of `AllAnchors`

For patterns with at most width w w w and depth d dd, the total runtime of `AllAnchors` is in

O(c w⋅d w 1/2).O\left(\frac{c^w \cdot d}{w^{1/2}}\right).O(w 1/2 c w⋅d​).

We restrict `Operations` on line 9 to only return the first d dd operations on the linear path in each direction, starting at the anchor operation: operations more than distance d dd away from the anchor cannot be part of a pattern of depth d dd.

We use the bound on the length of the list returned by calls to `AllAnchors` of Proposition 4.10 to bound the runtime. We can ignore the non-constant runtime of the concatenation of the outputs of recursive calls on line 35, as the total size of the outputs is asymptotically at worst of the same complexity as the runtime of the recursive calls themselves. Excluding the recursive calls, the only remaining lines of `AllAnchors` that are not executed in constant time are the `while` loop on lines 15–18 and the `Operations` and `sort` calls on lines 9–11. Using the same argument as in `CanonicalAnchors`, we can ignore the latter two calls by replacing the queue of operations by a lazy iterator of operations. The next operation given `op` and the graph `G` can always be computed in O(1)O(1)O(1) time using a depth-first traversal of `G`.

Consider the recursion tree of `AllAnchors`, i.e. the tree in which the nodes are the recursive calls to `AllAnchors` and the children are the executions spawned by the nested `for` loops on line 28–32. This tree has at most

C w=Θ(c w w 3/2)C_w = \Theta\left(\frac{c^w}{w^{3/2}}\right)C w​=Θ(w 3/2 c w​)

leaves. A path from the root to a leaf corresponds to a stack of recursive calls to `AllAnchors`. Along this recursion path, `seen_paths` set is always strictly growing (line 35) and the operations removed from `G` on lines 14 and 18 are all distinct. For each linear path, at most 2 d 2d2d operations are traversed. Thus the total runtime of the `while` loop (lines 15–18) along a path from root to leaf in the recursion tree is in O(w⋅d)O(w \cdot d)O(w⋅d). We can thus bound the overall complexity of executing the entire recursion tree by O(C w⋅w⋅d)=O(c w⋅d w 1/2)O(C_w \cdot w \cdot d) = O(\frac{c^w \cdot d}{w^{1/2}})O(C w​⋅w⋅d)=O(w 1/2 c w⋅d​).

1. Every anchor operation is on at least one previously unseen linear path, thus there can be at most w w w operations in the set of anchors. ↩︎

## 4.6. An automaton for multi-pattern matching

We have shown in the previous sections that graph pattern-matching can be reduced to a problem of tree inclusions, with trees of fixed width w w w. To complete the pattern-matching algorithm, we must provide a fast way to evaluate the subtree relation for many trees representing the set of all patterns we wish to match.

More precisely, for patterns P 1,…,P ℓ P_1, \dots, P_\ell P 1​,…,P ℓ​ with width w w w, fix a root operation r i r_i r i​ in P i P_i P i​ for each 1⩽i⩽ℓ 1 \leqslant i \leqslant \ell1⩽i⩽ℓ and consider the rooted tree duals of the canonical PSGs τ r i(P i π i)\tau_{r_i}(P_i^{\pi_i})τ r i​​(P i π i​​), with π i=π r i(P i)\pi_i = \pi_{r_i}(P_i)π i​=π r i​​(P i​) the canonical anchors. Then given a subject graph G GG, we wish to compute the set

{1⩽i⩽ℓ∣τ r i(P i π i)⊆τ r(G π)},\{1 \leqslant i \leqslant \ell \mid \tau_{r_i}(P_i^{\pi_i}) \subseteq \tau_r(G^\pi)\},{1⩽i⩽ℓ∣τ r i​​(P i π i​​)⊆τ r​(G π)},

for all anchor sets π∈Π r w(G)\pi \in \Pi_r^w(G)π∈Π r w​(G) and root operation r r r in G GG. This corresponds to the `IsSubTree` predicate introduced in the sketch of the algorith in section 4.4.

Instead of considering the trees of PSGs, it will prove easier to consider the contracted PSGs (cPSGs)

τ r i(c(P i π i))and τ r(c(G π)).\tau_{r_i}(c(P_i^{\pi_i}))\quad\textrm{and}\quad \tau_r(c(G^\pi)).τ r i​​(c(P i π i​​))and τ r​(c(G π)).

Such tree inclusions are equivalent to finding embeddings in the subject graph itself, provided that we keep track of the s p l i t splits pl it and c o n t r a c t contractcont r act weight maps (see section 4.3).

It will be useful to remind ourselves the following properties of contracted PSGs. Every operation of a cPSG (and thus every node in its rooted dual tree) is an anchor operation of the PSG. Per Proposition 4.4, the rooted dual tree of a cPSG is a ternary tree and has exactly w i d t h(G)−1 width(G) - 1 w idth(G)−1 nodes. Finally, recall the concept of an open value of a graph, i.e. a value that is missing either a use or define operation (see section 4.2).

#### Reduction of tree inclusion to string prefix matching #

Now consider two contracted spanning tree reductions c(G 1 π 1)c(G_1^{\pi_1})c(G 1 π 1​​) and c(G 2 π 2)c(G_2^{\pi_2})c(G 2 π 2​​) with values V 1 V_1 V 1​ and V 2 V_2 V 2​. To simplify notation, define

τ 1=τ r 1(c(G 1 π 1))and τ 2=τ r 2(c(G 2 π 2))\tau_1 = \tau_{r_1}(c(G_1^{\pi_1})) \quad\textrm{and}\quad \tau_2 = \tau_{r_2}(c(G_2^{\pi_2}))τ 1​=τ r 1​​(c(G 1 π 1​​))and τ 2​=τ r 2​​(c(G 2 π 2​​))

for some choice of root operations r 1 r_1 r 1​ and r 2 r_2 r 2​ in G 1 G_1G 1​ and G 2 G_2G 2​, respectively. We lift the ⊆\subseteq⊆ relation on rooted dual trees of PSGs introduced in section 4.5 to rooted dual trees of cPSGs in Such a way that there is an inclusion relation between two rooted dual trees of PSGs if and only if the same relation holds on the rooted duals of cPSGs.

We say that τ 1⊆τ 2\tau_1 \subseteq \tau_2 τ 1​⊆τ 2​ if and only if

- the trees share the same root operation,
- τ 1\tau_1 τ 1​ is a subtree of τ 2\tau_2 τ 2​,
- the s p i l t spiltspilt map coincides on the common subtree, and
- the c o n t r a c t contractcont r act map satisfies for all v∈V 1 v \in V_1 v∈V 1​:
  {c o n t r a c t(v)⊆c o n t r a c t(f(v))if v is an open value, c o n t r a c t(v)=c o n t r a c t(f(v))otherwise,\begin{cases}contract(v) \subseteq contract(f(v))\quad&\textrm{if }v\textrm{ is an open value},\\contract(v) = contract(f(v))\quad&\textrm{otherwise},\\\end{cases}{cont r act(v)⊆cont r act(f(v))cont r act(v)=cont r act(f(v))​if v is an open value, otherwise,​

where f: V 1↪V 2 f: V_1 \hookrightarrow V_2 f: V 1​↪V 2​ designates the embedding of V 1 V_1 V 1​ into V 2 V_2 V 2​ given by the tree embedding.

The first three conditions are taken as-is from the ⊆\subseteq⊆ relation on non-contracted trees, whilst the fourth condition on the c o n t r a c t contractcont r act map is specific to contracted trees.

Using Proposition 4.2, there are at most 2 open values for each linear path in the graph, and thus at most 2⋅w 2 \cdot w2⋅w open values in a rooted dual tree of a cPSG of width w w w. For each such contracted rooted dual, we can thus define a contracted string tuple S=(s 1,…,s 2 w)∈(O∗)2 w S = (s_1, \dots, s_{2w}) \in (O^\ast)^{2w}S=(s 1​,…,s 2 w​)∈(O∗)2 w given by the values of the c o n t r a c t contractcont r act map evaluated in the (up to) 2 w 2w2 w open values1.

If c o n t r a c t∣C contract|_Ccont r act∣C​ is the restriction of c o n t r a c t contractcont r act to the domain of definition of non-open values of a cPSG, the fourth condition for the inclusion relation ⊆\subseteq⊆ on rooted dual cPSGs, given above becomes an equality condition when restricted to non-open values. A special case of this property of particular interest to us is stated as the following result. The ⊆\subseteq⊆ relation on strings refers to prefix inclusion, i.e. s⊆t s \subseteq ts⊆t if and only if s ss is a prefix of t tt.

Proposition 4.14Inclusion of equal-width trees

Let S=(s 1,…,s 2 w)S = (s_1, \dots, s_{2w})S=(s 1​,…,s 2 w​) and T=(t 1,…,t 2 w)∈(O∗)2 w T = (t_1, \dots, t_{2w}) \in (O^\ast)^{2w}T=(t 1​,…,t 2 w​)∈(O∗)2 w be the contracted string tuples of τ 1\tau_1 τ 1​ and τ 2\tau_2 τ 2​ respectively. Then τ 1⊆τ 2\tau_1 \subseteq \tau_2 τ 1​⊆τ 2​ if and only if the trees share the same root, are isomorphic, have the same s p l i t splits pl it and c o n t r a c t∣C contract|_Ccont r act∣C​ maps and for all i∈{1,…,2 w}i \in \{1, \dots, 2w\}i∈{1,…,2 w}: s i⊆t i s_i \subseteq t_is i​⊆t i​.

The proof of this follows directly from observing that rooted duals of cPSGs have the same set of nodes and that the restriction to non-open values c o n t r a c t∣C contract|_Ccont r act∣C​ must satisfy equality.

Why restricting ourselves to trees of the same width w w w? It is sufficient for our purposes! All patterns are of width w w w by assumption and so are the rooted dual trees of the form τ r(G π)\tau_r(G^\pi)τ r​(G π), given that π∈Π r w(G)\pi \in \Pi_r^w(G)π∈Π r w​(G).

The string prefix matching problem is a simple computational task that can be generalised to check for multiple string patterns at the same time using a prefix tree. An overview of this problem can be found in appendix A. We can thus obtain a solution for the pattern matching problem for ℓ\ellℓ patterns:

Proposition 4.15Fixed anchor pattern matching

As above, let

- G GG be a graph, with π∈Π r w(G)\pi \in \Pi_r^w(G)π∈Π r w​(G) a set of w−1 w - 1 w−1 operations and r∈π r \in \pi r∈π a choice of root operation,
- P 1,…,P ℓ P_1, \dots, P_\ell P 1​,…,P ℓ​ be patterns of width w w w and depth d dd, with choices of root operations r 1,…,r ℓ r_1, \dots, r_\ell r 1​,…,r ℓ​ and canonical anchors π i=π r i(P i).\pi_i = \pi_{r_i}(P_i).π i​=π r i​​(P i​).

The set of all pattern embeddings mapping the canonical anchor set π i\pi_i π i​ to π\pi π and root r i r_i r i​ to r r r for 1≤i≤ℓ 1 \leq i \leq \ell1≤i≤ℓ can be computed in time O(w⋅d)O(w\cdot d)O(w⋅d) using at most ℓ\ellℓ pre-computed prefix tree of size at most (ℓ⋅d)w+1(\ell \cdot d)^w + 1(ℓ⋅d)w+1, each constructed in time complexity O((ℓ⋅d)w)O((\ell \cdot d)^w)O((ℓ⋅d)w).

For each pattern, we consider its canonical spanning tree reduction and construct a multi-dimensional prefix tree (see Appendix ) for each group of patterns that share the same spanning tree reduction.

Given a graph G GG, we can compute the cPSG of G GG for anchors π\pi π and map its rooted dual tree to the corresponding prefix tree. This can be done in O(∣T G∣)O(|T_G|)O(∣T G​∣) time by using a search tree. We can restrict G π G^\piG π to a graph of size O(w⋅d)O(w \cdot d)O(w⋅d) by truncating the linear paths to at most 2 d 2d2d length, as in the proof of Proposition 4.12. Thus we can assume ∣G π∣∈O(w⋅d)|G^\pi| \in O(w \cdot d)∣G π∣∈O(w⋅d).

The rest of the proof and the runtime follow from the multi-dimensional prefix tree construction detailed in Appendix ).

#### Combining everything #

Finally, putting Proposition . and Proposition 4.12 together, we obtain our main result.

Proposition 4.13Pattern matching

Let P 1,…,P ℓ P_1, \dots, P_\ell P 1​,…,P ℓ​ be patterns with width w w w and depth d dd. The pre-computation runs in time and space complexity

O((d⋅ℓ)w⋅ℓ+ℓ⋅w⋅d).O \left( (d\cdot \ell)^w \cdot \ell + \ell \cdot w \cdot d \right).O((d⋅ℓ)w⋅ℓ+ℓ⋅w⋅d).

For any subject graph G GG, the pre-computed prefix tree can be used to find all pattern embeddings P i→G P_i \to G P i​→G in time

O(∣G∣⋅c w w 1/2⋅d)O \left( |G| \cdot \frac{c^w}{w^{1/2}} \cdot d \right)O(∣G∣⋅w 1/2 c w​⋅d)

where c=6.75 c = 6.75c=6.75 is a constant.

The pre-computation consists of running the `CanonicalAnchors` procedure on each of the ℓ\ellℓ patterns and then transforming them into a map of prefix trees using Proposition .. By Proposition 4.7, `CanonicalAnchors` runs in O(w⋅d)O(w\cdot d)O(w⋅d) for each pattern, where we used that ∣P i∣⩽w⋅d|P_i| \leqslant w \cdot d∣P i​∣⩽w⋅d for all patterns. The total runtime of prefix construction is thus

O((d⋅ℓ)w⋅ℓ+ℓ⋅w⋅d).O \left( (d\cdot \ell)^w \cdot \ell + \ell \cdot w \cdot d \right).O((d⋅ℓ)w⋅ℓ+ℓ⋅w⋅d).

The complexity of pattern matching itself on the other hand is composed of two parts: the computation of all possible anchor sets Π r w(G)\Pi_r^w(G)Π r w​(G), and the execution of the prefix string matcher for each of the trees resulting from these sets π∈Π r w(G)\pi \in \Pi_r^w(G)π∈Π r w​(G). As `AllAnchors` must be run for every choice of root vertex r r r in G GG, the runtime is thus obtained by multiplying i) ∣G∣|G|∣G∣ with ii) the runtime of the prefix tree matching (Proposition .), and with iii) ∣Π r w(G)∣|\Pi_r^w(G)|∣Π r w​(G)∣, i.e. the number of anchor lists returned by `AllAnchors` (Proposition 4.10):

O(∣G∣⋅w⋅d⋅C w),O(|G| \cdot w \cdot d \cdot C_w ),O(∣G∣⋅w⋅d⋅C w​),

where C w C_w C w​ is the bound for the number of anchor lists returned by `AllAnchors`. The result follows.

1. The values can be ordered as usual by using the total lexicographic order on port labels of the tree. ↩︎

## 4.7. Benchmarks

Proposition 4.13 shows that pattern-independent matching can scale to large datasets of patterns but imposes some restrictions on the patterns and embeddings that can be matched. In this section, we discuss these limitations and give empirical evidence that the pattern-matching approach we have presented can be used on a large scale and outperform existing solutions.

#### Pattern limitations #

In section 4.2, we imposed conditions on the pattern embeddings to obtain a complexity bound for pattern-independent matching. We argued how these restrictions are natural for applications in quantum computing, and most of the arguments will also hold for a much broader class of computation graphs.

In future work, it would nonetheless be of theoretical interest to explore the importance of these assumptions and their impact on the complexity of the problem. As a first step towards a generalisation, our implementation and all our benchmarks in this section do not make any of these simplifying assumptions. Our results below give empirical evidence that a significant performance advantage can be obtained regardless.

#### Implementation #

We provide an open-source implementation in Rust of pattern independent matching using the results of this chapter, available on GitHub. The code and datasets used for the benchmarks themselves are available in a dedicated repository.

The implementation works for weighted or unweighted port graphs – of which typed minIR graphs are a special case – and makes none of the simplifying assumptions employed in the theoretical analysis. Pattern matching proceeds in two phases: precomputation and runtime.

Precomputation.In a first step, all graph patterns are processed and compiled into a single state automaton that will be used at runtime for fast pattern independent matching. The automaton in the implementation combines in one data structure two distinct computations of this chapter:

- the recursive branching logic used in the `AllAnchors` procedure to enumerate all possible choices of anchors.
- the automaton described in section 4.6 that matches patterns for a fixed set of anchors, and

The former is implemented with non-deterministic state transitions – each transition corresponding to choosing an additional anchor – , whereas the latter is implemented deterministically.

Concretely, the automaton is constructed by following the construction of section 4.4 to decompose each pattern into its canonical path-split graph. We then order the nodes of the PSG and express each node as a condition that ensures the connectivity and node weight in the graph matches the pattern. We thus obtain a chain of conditions, with a transition between any two consecutive conditions; transitions are deterministic by default and marked as non-deterministic whenever they lead to a condition on an anchor node. The state automaton for all patterns is then obtained by joining all chains of conditions into a tree.

Runtime.Pattern matching is then as simple as simulating the state automaton, evaluating all conditions on the graph G GG passed as input. The states in the automaton corresponding to the last condition of a pattern must be marked as end states, along with a label identifying the pattern that was matched. This can then be used at runtime to report all patterns found.

Our implementation has been tested for correctness, i.e. on the one hand that all matches that are reported are correct, and on the one hand that all pattern matches are found. This was done by comparing the matches of our implementation with the results obtained from matching every pattern separately on millions of randomly generated graphs and edge cases. We also ensured during benchmarking that the number of matches reported by our implementation and by Quartz were always the same.

#### Benchmarks #

Baseline.To assess practical use, we have benchmarked our implementation against a leading C++ implementation of pattern matching for quantum circuits from the Quartz superoptimiser project Xu, 2022Mingkuan Xu, Zikun Li, Oded Padon, Sina Lin, Jessica Pointing, Auguste Hirth, Henry Ma, Jens Palsberg, Alex Aiken, Umut A. Acar and Zhihao Jia. 2022. Quartz: Superoptimization of Quantum Circuits. In Proceedings of the 43rd ACM SIGPLAN International Conference on Programming Language Design and Implementation, June 2022. Association for Computing Machinery, 625--640. doi: 10.1145/3519939.3523433. This implementation is the principal component of an end-to-end quantum circuit optimisation pipeline. The results and speedups we obtain here thus apply and transfer directly to this application.

Dataset.We further ensure that our results apply in practice by using a real-world dataset of patterns. The Quartz optimiser finds opportunities for circuit optimisation by relying on precomputed equivalence classes of circuits (ECC). These are obtained exhaustively by enumerating all possible small quantum circuits, computing their unitaries and clustering them into classes of circuits with identical unitaries.

The generation of ECC sets is parametrised on the number of qubits, the maximum number of gates and the gate set in use. For these benchmarks we chose the minimal set of gates T, H, C X T, H, CX T, H, CX and considered circuits with up to 6 gates and 2, 3 or 4 qubits. The size of these pattern circuits is typical for the application1.

Thus, for our patterns, we have the bound d≤6 d \leq 6d≤6 for the maximum depth and width w=2, 3, 4 w = 2,3,4 w=2, 3, 4. In all experiments, the graph G GG subject to pattern matching was `barenco_tof_10` input, i.e. a 19-qubit circuit input with 674 gates obtained by decomposing a 10-qubit Toffoli gate using the Barenco decomposition Barenco, 1995Adriano Barenco, Charles H. Bennett, Richard Cleve, David P. DiVincenzo, Norman Margolus, Peter Shor, Tycho Sleator, John A. Smolin and Harald Weinfurter. 1995. Elementary gates for quantum computation. Physical Review A 52, 5 (November 1995, 3457--3467). doi: 10.1103/PhysRevA.52.3457.

Results.We study the runtime of our implementation as a function of the number ℓ\ellℓ of patterns being matched, up to ℓ=1 0 4\ell = 10^4ℓ=10 4 patterns. We expect the runtime of pattern matching algorithms that match one pattern at a time to scale linearly with ℓ\ellℓ. On the other hand, Proposition 4.13 results in a complexity that is independent of ℓ\ellℓ.

For each value of ℓ\ellℓ, we select a subset of all patterns in the ECC sets at random. For w=2 w = 2 w=2, there are only a total of ℓ=1954\ell = 1954ℓ=1954 patterns, explaining why we do not report result beyond that number. For ℓ=200\ell = 200ℓ=200 patterns, our proposed algorithm is 3×3\times3× faster than Quartz. As expected, the advantage of our approach increases as we match more patterns, scaling up to a 20×20\times20× speedup for ℓ=1 0 4\ell=10^4ℓ=10 4. The results are summarised in the following figure:

Runtime of pattern matching for ℓ=0…1 0 4\ell = 0\dots 10^4ℓ=0…10 4 patterns on 2, 3 and 4 qubit quantum circuits from the Quartz ECC dataset, for our implementation (Portmatching) and the Quartz project. All ℓ=1954\ell = 1954ℓ=1954 two-qubit circuits were used, whereas for 3 and 4 qubit circuits, ℓ=1 0 4\ell = 10^4ℓ=10 4 random samples were drawn.

Dependency on w w w and ℓ\ellℓ.We further study the runtime of our algorithm as a function of its two main parameters, the number of patterns ℓ\ellℓ and the pattern width w w w, on an expanded dataset. To this end, we generate random sets of 10,000 pattern circuits with 15 gates and between w=2 w=2 w=2 and w=10 w=10 w=10 qubits, using the same gate set as previously. The resulting pattern matching runtimes are shown in the figure below.

From Proposition 4.13, we expect that the pattern matching runtime is upper bounded by a ℓ\ellℓ-independent constant. Our results support this result for w=2 w=2 w=2 and w=3 w=3 w=3 qubit patterns, where runtime seems indeed to saturate, reaching an observable runtime plateau at large ℓ\ellℓ.

We suspect on the other hand that the exponential c w c^wc w dependency in the complexity bound of Proposition 4.13 makes it difficult to observe a similar plateau for w≥4 w \geq 4 w≥4, as we expect this upper bound on the runtime to increase rapidly with qubit counts . A runtime ceiling is not directly observable at this experiment size, but the gradual decrease in the slope of the curve is consistent with the existence of the ℓ\ellℓ-independent upper bound predicted in Proposition 4.13.

Runtime of our pattern matching for random quantum circuits with up to 10 qubits.

1. Such small circuit sizes are imposed in part by the fact that ECCs of larger circuits quickly become unfeasible to generate as their number grows exponentially. In practice, large circuit transformations can often be expressed as the composition of smaller atomic transformations, hence the good performance of this approach in practice. ↩︎