arXiv CSopen access

gMatch: Fine-Grained and Hardware-Efficient Subgraph Matching on GPUs

arXiv:2604.10601v1 [cs.DB] 12 Apr 2026

gMatch: Fine-Grained and Hardware-Efficient Subgraph Matching on GPUs Weitian Chen

Shixuan Sun

Cheng Chen

Shanghai Jiao Tong University [email protected]

Shanghai Jiao Tong University [email protected]

ByteDance Inc [email protected]

Yongmin Hu

Yingqian Hu

Minyi Guo

ByteDance Inc [email protected]

ByteDance Inc [email protected]

Guizhou University [email protected]

ABSTRACT Subgraph matching is a core operation in graph analytics, supporting a broad spectrum of applications from social network analysis to bioinformatics. Recent GPU-based approaches accelerate subgraph matching by leveraging parallelism but rely on a coarse-grained execution model that suffers from scalability and efficiency issues due to high memory overhead and thread underutilization. In this paper, we propose gMatch, a hardware-efficient subgraph matching approach on GPUs. gMatch introduces a fine-grained execution model that reduces memory consumption and enables flexible task scheduling among threads. We further design warp-level batch exploration and lightweight load balancing to improve execution efficiency and scalability. Experiments on diverse workloads and real-world datasets show that gMatch outperforms state-of-the-art subgraph matching methods, including STMatch, T-DFS, and EGSM, in both performance and scalability. We also compare against stateof-the-art systems for mining small patterns, such as BEEP and G2 Miner. While these systems achieve better performance on small datasets, gMatch scales to substantially larger queries and datasets, where existing approaches degrade or fail to complete. PVLDB Reference Format: Weitian Chen, Shixuan Sun, Cheng Chen, Yongmin Hu, Yingqian Hu, and Minyi Guo. gMatch: Fine-Grained and Hardware-Efficient Subgraph Matching on GPUs. PVLDB, 19(8): XXX-XXX, 2026. doi:XX.XX/XXX.XX PVLDB Artifact Availability: The source code, data, and/or other artifacts have been made available at https://github.com/SJTU-Liquid/gMatch.

1

INTRODUCTION

Given a query graph 𝑄 and a data graph 𝐺, subgraph matching aims to find all embeddings of 𝑄 in 𝐺. For example, in Figure 1, the mapping 𝑀 = {(𝑢 1, 𝑣 58 ), (𝑢 2, 𝑣 57 ), (𝑢 3, 𝑣 1 ), (𝑢 4, 𝑣 60 )} is a valid embedding, or match, of 𝑄 in 𝐺. As a fundamental operation in graph analysis, subgraph matching underpins a wide range of real-world This work is licensed under the Creative Commons BY-NC-ND 4.0 International License. Visit https://creativecommons.org/licenses/by-nc-nd/4.0/ to view a copy of this license. For any use beyond those covered by this license, obtain permission by emailing [email protected]. Copyright is held by the owner/author(s). Publication rights licensed to the VLDB Endowment. Proceedings of the VLDB Endowment, Vol. 19, No. 8 ISSN 2150-8097. doi:XX.XX/XXX.XX

𝑢3

𝑣1

𝑣2

𝑣50

C

C

… C

𝑣56

C

… C

C

𝑢2

C

B

B 𝑣57

A 𝑢1

𝑣51 𝑣52

D 𝑢4

(a) Query graph 𝑄.

A

B

𝑣58

𝑣59

D

D

D

𝑣60

𝑣61

𝑣62

(b) Data graph 𝐺.

Figure 1: Example query graph and data graph. applications, including social network recommendation [15], fraud detection [21], bioinformatics [6, 46], and cheminformatics [39]. Given the importance of subgraph matching, numerous algorithms have been proposed [5, 16, 33]. Although they adopt different optimization techniques, most follow the same core idea: iteratively extend partial matches by mapping query vertices in 𝑄 to data vertices in 𝐺 according to a predefined order. However, since subgraph isomorphism is NP-complete, the search space can be prohibitively large. To accelerate the computation, recent studies leverage GPUs to parallelize the matching process [34, 36, 38, 40, 43, 45]. All these methods adopt a coarse-grained parallel execution model for subgraph matching, where each partial match 𝑀 is treated as a parallel task and a warp, the basic scheduling unit on GPUs, serves as the worker. The task involves extending 𝑀 by mapping the next query vertex to valid data vertices. Based on their search strategies, these methods can be classified into two categories: breadth-first search (BFS)-based and depth-first search (DFS)-based. Early approaches such as GpSM [36], GSI [45], and cuTS [40] adopt a BFS-based search strategy, which explores the search space level by level. At each iteration, all partial matches are extended by mapping the next query vertex to candidate data vertices. This strategy enables a large number of parallel tasks, allowing efficient GPU utilization and balanced workload distribution across warps. However, it requires storing all partial matches at every level, leading to high memory consumption due to the large search space. As a result, the scalability of BFS-based methods is severely limited by GPU memory capacity. To reduce memory overhead, recent methods such as T-DFS [43], STMatch [38], and EGSM [34] adopt a DFS-based search strategy. Instead of expanding all partial matches level by level, each warp recursively explores the search space. Specifically, after generating new partial matches from a given match 𝑀, the warp selects one

to continue the search. By avoiding the need to store all partial matches at each level, this strategy significantly improves scalability compared to BFS-based approaches. Key Insights. Despite its widespread adoption, the coarse-grained parallel execution model fundamentally conflicts with real-world power-law graphs. Because each partial match has highly dynamic and skewed workloads, while GPUs adopt the SIMT execution model and a warp provides fixed and inflexible computing resources, this design introduces fundamental limitations in both scalability and computational efficiency for GPU-based subgraph matching. First, the execution stack incurs high memory overhead. To handle dynamic workloads, the coarse-grained method must pre-allocate a large execution stack of size 𝑂 (|𝑉 (𝑄)| ×𝑑 max ) for each warp, where 𝑑 max is the maximum degree of 𝐺. The stack stores intermediate results during the search, e.g., feasible candidates and derived partial matches. Since 𝑑 max can reach millions in power-law graphs and GPUs require hundreds of active warps for full utilization, this high memory cost limits scalability and reduces concurrency. Second, the mismatch between the fixed warp size (32 threads) and the dynamic workloads leads to severe thread underutilization. To extend a partial match, the coarse-grained method performs set intersections over the neighbors of previously matched data vertices. Most candidate vertices have very small degrees after filtering. A warp must still execute 32 threads in lock-step, causing pervasive idle threads and wasted GPU parallelism. These limitations are not caused by specific implementations; they stem from the inherent structure of the coarse-grained paradigm itself. Overcoming them requires fundamentally rethinking how subgraph matching should be parallelized on GPUs. Our Work. In this paper, we propose gMatch, a hardware-efficient subgraph matching approach on GPUs. To resolve the structural mismatch between the GPU execution model and the coarse-grained parallel execution, we first introduce a fine-grained parallel execution model that redefines the basic task and the basic worker. Instead of letting a warp process one partial match, we let a single thread process the lightweight task of checking whether a query vertex can be mapped to a data vertex. This design provides two key advantages. First, it reduces the per-worker state from 𝑂 (|𝑉 (𝑄)| × 𝑑 max ) to 𝑂 (|𝑉 (𝑄)|), eliminating stack-induced memory bottlenecks and enabling the system to process much larger graphs while supporting more concurrent warps. Second, it decomposes partial-match processing into naturally divisible units, allowing the workload of a single partial match to be distributed across multiple warps or tasks from different partial matches to be combined within a warp. However, scheduling such lightweight tasks at scale is challenging, because tens of thousands of threads execute concurrently and GPU threads are organized hierarchically into warps and thread blocks. To address this, we propose a warp-level batch exploration mechanism. We merge task groups (i.e., the tasks involved in processing a partial match) into a virtual task pool and let threads within a warp pull tasks using two shared pointers that track task boundaries. This mechanism provides minimal memory overhead, high warp utilization, and efficient composability across tasks. Such batch processing is impossible under prior coarse-grained models because each partial match is an indivisible unit; processing multiple partial matches together would require storing large intermediate states, leading to prohibitive memory costs (e.g., STMatch).

Load balancing is essential for efficient subgraph matching on GPUs. Our key observation is that although the root of the search tree may offer limited parallelism, the number of partial matches expands exponentially as the search proceeds, making the workload inherently embarrassingly parallel. When the number of tasks far exceeds the number of workers, load balancing becomes naturally easy because task variability is automatically smoothed out across a large shared task pool. Building on this observation, we introduce a lightweight load balancing strategy that includes an initialization phase to generate a large pool of partial matches as initial tasks. This creates abundant initial parallelism and makes load balancing largely self-sustaining, minimizing the need for work stealing and reducing scheduling overhead. Consequently, work stealing becomes a supplement role in gMatch rather than the primary mechanism, as in prior methods. We evaluate gMatch on two representative workloads: (1) large queries on medium-sized graphs with multiple labels, and (2) small queries on large graphs with few labels. Experimental results show that gMatch significantly outperforms state-of-the-art methods, including STMatch [38], T-DFS [43], and EGSM [34], while offering better scalability. We also compare our approach with systems designed for mining small patterns, including BEEP [23] and G2 Miner [8]. Although these systems outperform our approach on small datasets, our method scales to substantially larger queries and datasets.

2

BACKGROUND

In this section, we introduce the background related to this paper.

2.1

Preliminaries

We focus on undirected, labeled graphs 𝑔 = (𝑉 , 𝐸), where 𝑉 is the set of vertices and 𝐸 ⊆ 𝑉 × 𝑉 is the set of edges. For a vertex 𝑢 ∈ 𝑉 , let 𝑁 (𝑢) denote its neighbor set and 𝑑 (𝑢) = |𝑁 (𝑢)| its degree. Each vertex is associated with a label via a labeling function 𝐿 : 𝑉 → Σ, where Σ is the label set. We denote the query graph by 𝑄 and the data graph by 𝐺. Vertices and edges in 𝑄 are referred to as query vertices and query edges, and those in 𝐺 as data vertices and data edges. Frequently used notations are summarized in Table 1. Definition 2.1 introduces subgraph isomorphism, which we refer to as a match. The goal of subgraph matching is to find all matches of 𝑄 in 𝐺. Definition 2.1. Given graphs 𝑔 and 𝑔′ , a subgraph isomorphism is an injective function 𝑀 : 𝑉 (𝑔) → 𝑉 (𝑔′ ) s.t. 1) ∀𝑣 ∈ 𝑉 (𝑔), 𝐿(𝑣) = 𝐿(𝑀 [𝑣]); 2) ∀𝑒 (𝑣, 𝑣 ′ ) ∈ 𝐸 (𝑔), 𝑒 (𝑀 [𝑣], 𝑀 [𝑣 ′ ]) ∈ 𝐸 (𝑔′ ). Given 𝑄 and 𝐺 in Figure 1, {(𝑢 1, 𝑣 58 ), (𝑢 2, 𝑣 57 ), (𝑢 3, 𝑣 1 ), (𝑢 4, 𝑣 60 )} is a match. We can observe that there exist 112 distinct matches of 𝑄 in 𝐺. The subgraph matching problem requires finding all such mappings from 𝑄 to 𝐺.

2.2

Coarse-Grained Parallel Execution

Although existing methods apply various optimizations, they share the same core principle: given a query graph 𝑄 and a data graph 𝐺, they iteratively extend intermediate results along a matching order 𝜙 to find all matches. A matching order 𝜙 is a permutation of the query vertices, where 𝜙 [𝑖] denotes the 𝑖-th vertex in the sequence

Table 1: Frequently used symbols and notations.

𝜙 𝑣

𝑢

Notations

Descriptions

𝑔, 𝑄 and 𝐺 𝑉 (𝑔), 𝐸 (𝑔) and Σ(𝑔) 𝑑 (𝑢), 𝐿(𝑢) and 𝑁 (𝑢) 𝑒 (𝑢, 𝑢 ′ ) 𝑑 max, 𝑑 avg 𝐹 (𝑢) 𝐶𝑀 𝐿 (𝑢) 𝐶𝑀 𝜙 𝜙 𝑁 + (𝑢)

graph, query graph and data graph vertex set, edge set and label set degree, label and neighbor set of 𝑢 undirected edge connecting 𝑢 and 𝑢 ′ maximum and average degree of the graph feasible candidate set of 𝑢 given partial match 𝑀 local candidate set of 𝑢 given partial match 𝑀 matching order backward neighbors of 𝑢 given 𝜙

𝐶

𝑣

𝑢

𝑣

𝑣

…

𝑣 …

𝑣

𝑣

𝑣 𝑣

Output

𝑣

𝑣 …

𝑢 𝑣 𝑣 𝑣 𝑣

← 𝑁 𝑀[𝑢 ] ∩ 𝑁 𝑀[𝑢 ] =𝑁 𝑣 ∩𝑁 𝑣

𝑣 …

𝑢

𝑢

𝑣 𝑣

𝑣

𝑣

𝑣

Figure 2: The search tree of coarse-grained parallel execution on the example graphs in Figure 1.

𝜙

(1 ⩽ 𝑖 ⩽ |𝑉 (𝑄)|). For each query vertex 𝑢, let 𝑁 + (𝑢) denote the 𝜙 set of neighbors that appear earlier than 𝑢 in 𝜙. We call 𝑁 + (𝑢) the backward neighbors of 𝑢. The matching order is required to Ó 𝜙 be connected, i.e., ∀𝑢 ∈ 𝑉 (𝑄) 𝑢 ≠ 𝜙 [1], 𝑁 + (𝑢) ≠ ∅. Algorithm 1 illustrates the parallel subgraph matching process. Based on their traversal strategies, existing methods can be categorized into BFSbased and DFS-based. BFS-Search. Lines 1–7 describe a BFS-based approach. It begins by generating an initial set M of partial matches by mapping the first query edge (based on the matching order 𝜙) to data edges that satisfy vertex label constraints. A partial match with 𝑖 query vertices corresponds to a subgraph isomorphism from the induced subgraph of the first 𝑖 vertices in 𝜙 to 𝐺. Each 𝑀 ∈ M is then iteratively extended by mapping the next query vertex 𝑢 to candidate data vertices. At each iteration, M is processed in parallel, with each warp handling one partial match. Lines 19–26 perform the extension by mapping 𝑢 to data vertices 𝑣 that satisfy label and connectivity 𝜙 constraints, i.e., 𝑣 is adjacent to the data vertices mapped to 𝑁 + (𝑢). 𝐹 The set of such candidates is denoted by 𝐶𝑀 (𝑢) (i.e., the feasible candidate sets), and the extension is executed in parallel across threads within a warp, leveraging intra-warp parallelism. Once all query vertices are matched, results are reported (Line 7). GpSM [36] and GSI [45] adopt this BFS-based strategy. DFS-Search. Storing all partial matches at each iteration incurs high memory overhead and limits scalability. To address this, recent methods such as STMatch [38], T-DFS [43], and EGSM [34] employ a depth-first strategy (Lines 8–11), where each warp recursively extends a single partial match by mapping one query vertex at a time (Lines 12–16). This reduces memory usage, as only the current path in the search tree is maintained. Abstraction. The exploration process builds a search tree where each node represents a partial match 𝑀, and its children M𝑀 are the new partial matches obtained by mapping the next query vertex 𝐹 (𝑢). The difference lies in the 𝑢 ∈ 𝜙 to feasible candidates in 𝐶𝑀 traversal strategy: BFS-based methods explore the tree in a breadthfirst manner, while DFS-based methods use depth-first traversal. Despite this difference, both approaches adopt a coarse-grained parallel model, where each warp is assigned a partial match as a task and extends it independently. Intra-warp parallelism is employed to compute feasible candidates for extending each partial match. Take 𝑄 and 𝐺 in Figure 1 as an example. The subgraph matching algorithm first generates a matching order 𝜙 = (𝑢 1, 𝑢 2, 𝑢 3, 𝑢 4 ) and then performs a BFS or DFS search following 𝜙, as shown in the

search tree (Figure 2). At level 𝑖 of the tree, the algorithm computes feasible candidates for query vertex 𝜙 [𝑖] based on the current partial match (i.e., the path from the root). The path (𝑣 58, 𝑣 59 ) corresponds to a match of the subgraph induced by {𝑢 1, 𝑢 2 }. To extend this match to 𝑢 3 , the algorithm computes the intersection 𝑁 (𝑣 58 ) ∩ 𝑁 (𝑣 59 ) and retains only vertices with the same label as 𝑢 3 . This yields 𝐹 (𝑢 ) = {𝑣 , 𝑣 , ..., 𝑣 }. Note that for the feasible candidate set 𝐶𝑀 3 51 52 56 GPU-based subgraph matching, the set intersection operation is computed by a warp, where each thread is assigned a vertex in 𝑁 (𝑣 59 ) and checks the existence of this vertex in 𝑁 (𝑣 58 ). Remark. Our paper, together with prior work [32, 33, 47], identifies two representative workload patterns: (1) large queries over medium-sized graphs, and (2) small queries over large-scale graphs. A detailed discussion on these workloads, including why existing coarse-grained methods struggle to generalize across both efficiently, is deferred to Appendix C due to space limitations.

2.3

Issues in Coarse-Grained Parallel Execution

Although DFS-based search significantly reduces memory consumption compared to BFS-based search, we observe that its coarsegrained parallel execution leads to severe performance bottlenecks. Issue 1: The execution stack incurs high memory overhead, limiting both scalability and computational efficiency. The execution stack stores intermediate results during the search process. As shown in Algorithm 1, when processing a partial match 𝑀, a 𝐹 (𝑢) and the set of warp must store both the feasible candidates 𝐶𝑀 𝐹 (𝑢) derived partial matches M𝑀 . Since the maximum size of 𝐶𝑀 and M𝑀 is bounded by the maximum degree 𝑑 max of 𝐺, the coarsegrained strategy requires a buffer of size 𝑂 (𝑑 max ). Dynamic memory 𝐹 (𝑢) and M is impractical allocation based on the actual size of 𝐶𝑀 𝑀 on GPUs due to high allocation overhead. Given that the DFS stack depth is |𝑉 (𝑄)|, the total space complexity of the execution stack per warp is 𝑂 (|𝑉 (𝑄)| × 𝑑 max ). Real-world graphs typically follow a power-law degree distribution, where a small fraction of vertices have a large number of neighbors. This leads to high memory overhead for the execution stack, limiting scalability and efficiency on GPUs, which offer abundant compute resources but limited memory. For example, in LDBC, a social network benchmark, 𝑑 max reaches 4 million, requiring 16 MB to store feasible candidates per partial match (assuming 4-byte vertex IDs). On an NVIDIA 4090 GPU with 128 streaming multiprocessors (SMs), each capable of 64 active

Global Memory

Procedure BFS-Search(𝑄, 𝐺, 𝜙) Ó 2 M ← { { (𝜙 [1], 𝑣), (𝜙 [2], 𝑣 ′ ) } |𝑒 (𝑣, 𝑣 ′ ) ∈ 𝐸 (𝐺 ) 𝐿 (𝜙 [1] ) = Ó ′ 𝐿 (𝑣) 𝐿 (𝜙 [2] ) = 𝐿 (𝑣 ) }; 3 for 𝑖 ← 3 to |𝜙 | do 4 M ′ ← M, M ← ∅; /* Inter-warp parallelism. */ 5 parallel foreach 𝑀 ∈ M ′ do Ð 6 M ← M Process(𝐺, 𝑀, 𝜙 [𝑖 ]);

Execution Stack of a Warp

1

7

Output M;

Buffers 𝑢1

𝑣58

𝑢2

𝑣57 𝑣59

𝑢3

𝑣1 𝑣2 𝑣3 … 𝑣50

𝑢4

18

Procedure Process(𝐺, 𝑀, 𝑢) 𝐹 (𝑢 ) ← { }; 𝐶𝑀 /* Intra-warp parallelism.

20 21 22

24 25 26

…

Warp #N-1

Idle Threads

Execution Stack of a Warp Buffers 𝑢1

𝑣58

𝑢2

𝑣57 𝑣59

𝑢3

𝑣51 𝑣52 … 𝑣56

Thread Block Warp #0

Warp #1

Partial Match 𝑣58 𝑣59

Warp #2

𝑁 𝑣59 𝑣51 𝑣52 … 𝑣56 𝑣58 𝑣60 𝑣61

𝑢4

…

…

Warp #N-1

Idle Threads

Figure 3: Coarse-grained parallel execution on the example graphs in Figure 1.

*/

𝐿 (𝑢 ) ← 𝑁 (𝑀 [𝑢 ′ ] ) where 𝑢 ′ ∈ 𝑁 (𝑢 ); 𝐶𝑀 + 𝐿 (𝑢 ) do parallel foreach 𝑣 ∈ 𝐶𝑀 Ó if 𝐿 (𝑢 ) = 𝐿 (𝑣) 𝑣 is not mapped in 𝑀 then 𝐹 ← 𝑡𝑟𝑢𝑒;

foreach 𝑢 ′′ ∈ 𝑁 + (𝑢 ) where 𝑢 ′′ ≠ 𝑢 ′ do if 𝑣 ∉ 𝑁 (𝑀 [𝑢 ′′ ] ) then 𝐹 ← 𝑓 𝑎𝑙𝑠𝑒, break; 𝐹 (𝑢 ) ← 𝐶 𝐹 (𝑢 ) Ð {𝑣 }; if 𝐹 = 𝑡𝑟𝑢𝑒 then 𝐶𝑀 𝑀 Ð 𝐹 return {𝑀 { (𝑢, 𝑣) } |𝑣 ∈ 𝐶𝑀 (𝑢 ) }; 𝜙

23

…

(b) Extending partial match { (𝑢 1 , 𝑣58 ), (𝑢 2 , 𝑣59 ) }.

𝜙

19

…

…

17

Warp #2

𝑁 𝑣57

(a) Extending partial match { (𝑢 1 , 𝑣58 ), (𝑢 2 , 𝑣57 ) }.

𝑑max

Procedure Search(𝑄, 𝐺, 𝑀, 𝜙, 𝑖) 13 if 𝑖 = |𝜙 | + 1 then Output 𝑀, return; 14 M𝑀 ← Process(𝐺, 𝑀, 𝜙 [𝑖 ]); 15 foreach 𝑀 ′ ∈ M𝑀 do 16 Search(𝑄, 𝐺, 𝑀 ′ , 𝜙, 𝑖 + 1);

Warp #1

Partial Match 𝑣58 𝑣57

𝑑max

8

12

Warp #0

𝑣1 𝑣2 𝑣3 … 𝑣50 𝑣58 𝑣60 𝑣61

Global Memory

Procedure DFS-Search(𝑄, 𝐺, 𝜙) Ó 9 M ← { { (𝜙 [1], 𝑣), (𝜙 [2], 𝑣 ′ ) } |𝑒 (𝑣, 𝑣 ′ ) ∈ 𝐸 (𝐺 ) 𝐿 (𝜙 [1] ) = Ó ′ 𝐿 (𝑣) 𝐿 (𝜙 [2] ) = 𝐿 (𝑣 ) }; /* Inter-warp parallelism. */ 10 parallel foreach 𝑀 ∈ M do 11 Search(𝑄, 𝐺, 𝑀, 𝜙, 3);

Thread Block

…

Algorithm 1: Coarse-Grained Parallel Subgraph Matching

warps and 24 GB memory, the execution stack quickly becomes a bottleneck. For a query graph with just five vertices and eight active warps per SM, the stack alone would exceed 80 GB of memory. Reducing the number of active warps to save memory significantly underutilizes computational resources. Additionally, the stack must reside in global memory (as it is too large to fit in shared memory), further degrading data access speed. Figure 3a illustrates the issue of high memory consumption of the execution stack, using the example graph of Figure 1. To extend the partial match 𝑀 = {(𝑢 1, 𝑣 58 ), (𝑢 2, 𝑣 57 )}, the coarse-grained parallel approach computes the intersection 𝑁 (𝑣 58 ) ∩ 𝑁 (𝑣 57 ) and 𝐹 (𝑢 ) = filters out vertices whose label is not 𝐿(𝑢 3 ). The result 𝐶𝑀 3 {𝑣 1, 𝑣 2, ..., 𝑣 50 } is completely computed and materialized in the third level of the stack. The result has 50 vertices and is stored in the preallocated buffer of size 𝑂 (𝑑 max ). This preallocated buffer incurs high memory overhead. Issue 2: The mismatch between the fixed warp size and the dynamic workload causes severe thread underutilization. The coarse-grained parallel execution assigns one partial match to each warp, where 32 threads execute in lockstep. However, as shown in

𝐿 (𝑢) Lines 17–26 of Algorithm 1, the size of the local candidate set 𝐶𝑀 ′ (i.e., 𝑁 (𝑀 [𝑢 ])) varies dynamically during the search. Real-world graphs typically follow a power-law degree distribution, where 𝐿 (𝑢)| is often most vertices have few neighbors. As a result, |𝐶𝑀 smaller than the warp size, leading to poor intra-warp parallelism and underutilized threads. Figure 3b further illustrates this issue, using the example graphs of Figure 1. When extending the partial match 𝑀 = {(𝑢 1, 𝑣 58 ), (𝑢 2, 𝑣 59 )}, the coarse-grained parallel approach assigns each thread 𝐿 (𝑢 ), namely 𝑁 (𝑣 ). Each lane within the warp a vertex from 𝐶𝑀 3 59 then check the validity of the vertex and thus construct the feasible candidates for the next level. However, there are only 9 vertices within 𝑁 (𝑣 59 ), and 23 threads within the warp are idle. Moreover, this underutilization problem worsens at the next level. This is because each feasible candidate identified at the current level must then verify its own local candidate set (which is also 𝑁 (𝑣 59 )), further amplifying the inefficiency. STMatch [38] addresses this issue using loop unrolling, which sets a fixed unrolling factor 𝜎 (e.g., 2, 4, or 8) to process multiple partial matches per warp by unrolling the loop at Line 15 in Algorithm 1. While loop unrolling reduces the idle rate, it does not fully 𝐿 (𝑢)| varies resolve thread underutilization since 𝜎 is fixed while |𝐶𝑀 dynamically. More critically, it increases memory consumption by expanding the execution stack to 𝑂 (|𝑉 (𝑄)| × 𝑑 max × 𝜎), further limiting scalability. To quantify this thread underutilization issue, we define the 𝐿 (𝑢 ) | 32− |𝐶𝑀 1 Í idle rate as | M , where M is the set of partial 𝑀∈M 32 | matches generated during the search. We generate 100 random 12-vertex query graphs per dataset and obtain the average idle rate. Table 2 shows the idle rate across five datasets, revealing severe hardware underutilization (see Table 3 for details of these datasets).

Table 2: Idle rates with varying loop unrolling sizes. Unrolling Size

db

en

gw

gh

wt

1 2 4 8

70.74% 50.22% 30.16% 19.32%

45.14% 31.38% 20.75% 10.84%

40.34% 38.62% 22.94% 15.54%

48.19% 59.47% 21.62% 11.41%

58.80% 29.82% 12.49% 9.80%

Summary. Existing methods suffer from significant performance limitations due to their coarse-grained parallel execution model on GPUs: (1) high memory overhead and poor scalability, as the total space complexity required for execution stacks grows as 𝑂 (|𝑉 (𝑄)|× 𝑑 max ); and (2) severe thread underutilization, as the neighbor sets explored by a warp are often much smaller than the warp size, leaving many threads idle.

3

AN OVERVIEW OF GMATCH Input Query Graph 𝑄

Data Graph 𝐺

Thread Block Warp #0

Lane 𝐿 𝑣1′ 𝐶𝑀 ′

#30 #31

Warp #1 Warp #2

…

… … … …

𝑣𝑛′′

…

𝑣1′′ 𝑣2′′ 𝐿 𝐶𝑀 ′′

…

…

Initial Tasks

…

′ 𝑣𝑚

Global Memory Data Graph

…

Partial Matches

𝑀′ 𝑀′′

Matching Order 𝜙

#0 #1

𝑣2′

Warp #N-1

Shared Memory Warp #0 Warp #1

…

…

Auxiliary Data Structure

Warp #N-1

Local Task Queue

Figure 4: An overview of gMatch. Figure 4 presents an overview of gMatch, a hardware-efficient subgraph matching approach. Given 𝑄 and 𝐺, we first generate a matching order 𝜙. As matching order selection has been extensively studied and its overhead is negligible, we generate 𝜙 on the CPU using the RI method [6], which has demonstrated strong performance in prior work [32, 47]. gMatch is compatible with other matching order generation strategies as well. gMatch focuses on efficient parallel execution on GPUs. Unlike prior coarse-grained parallel model that treat a warp as the basic unit of work, we propose a fine-grained parallel execution model that treats each thread as a worker and each candidate validation as a task. This model enables fine-grained processing of partial matches and performs the search on-the-fly without materializing intermediate results, significantly reducing execution stack space and eliminating scalability constraints. Building on this model, we introduce a warp-level batch exploration technique that groups multiple partial matches into a single warp. This aligns with the GPU architecture, where a warp is the basic scheduling unit. The technique eliminates thread underutilization without incurring additional memory overhead. To further improve scalability, we design a lightweight load balancing strategy based on the observation that

subgraph matching is an embarrassingly parallel problem. At the start of execution, we perform a BFS-based exploration to generate partial matches until a specified threshold is reached, forming an initial task pool. Warps fetch tasks independently from this pool. Once the pool is exhausted, idle warps perform work stealing to maintain balance. Thanks to the large initial task pool, stealing is rarely triggered, keeping scheduling overhead low. We detail each of these components in the following sections.

4

EFFICIENT PARALLEL EXECUTION

In this section, we present the efficient parallel execution strategy employed in gMatch.

4.1

Fine-Grained Parallel Execution

Given 𝑄, 𝐺, and 𝜙, existing methods [34, 38, 40, 43] treat each partial match 𝑀 as a parallel task and assign it to a warp. Suppose that 𝑀 contains 𝑖 matched vertices and 𝑢 is the (𝑖 + 1)-th query vertex in 𝜙. The warp extends 𝑀 by mapping 𝑢 to candidates 𝑣 in the feasible 𝐹 (𝑢). To store 𝐶 𝐹 (𝑢) and the resulting partial matches, each set 𝐶𝑀 𝑀 warp allocates a buffer of size 𝑂 (𝑑 max ). As a result, a worker must maintain an execution stack of size 𝑂 (𝑑 max × |𝑉 (𝑄)|), leading to significant memory overhead as discussed in Section 2.3. To address these limitations, we design a fine-grained execution model. Instead of assigning a partial match 𝑀 to a warp as a monolithic task, we decompose it into smaller and independent tasks of the form 𝑇𝑀 (𝑢, 𝑣), representing an attempt to extend 𝑀 by mapping 𝑢 to a candidate vertex 𝑣. A task is valid if 𝑣 ∈ 𝑁 (𝑀 [𝑢 ′ ]) 𝜙 for all 𝑢 ′ ∈ 𝑁 + (𝑢). To reduce candidate enumeration, we fur𝐿 (𝑢) = 𝑁 (𝑀 [𝑢 ′ ]), ther restrict attention to a local candidate set 𝐶𝑀 where 𝑢 ′ = arg min𝑢 ′ ∈𝑁 𝜙 (𝑢 ) |𝑁 (𝑀 [𝑢 ′ ])|, i.e., the neighbor set of the + mapped backward query vertex with the fewest neighbors. Thus, we define the corresponding set of tasks as the task group for 𝑀, 𝐿 (𝑢)}. denoted by 𝑇𝑀𝐺 (𝑢) = {𝑇𝑀 (𝑢, 𝑣) | 𝑣 ∈ 𝐶𝑀 In GPUs, threads are the fundamental execution units, while warps serve as the basic scheduling units. Our execution model assigns each task group 𝑇𝑀𝐺 (𝑢) to a warp, and distributes its tasks 𝑇𝑀 (𝑢, 𝑣) to individual threads within the warp. Each thread independently processes its assigned task by verifying the feasibility of mapping 𝑢 to 𝑣. After task execution, the warp gathers the valid candidates and generates the corresponding extended partial matches. Figure 5 illustrates the search tree in fine-grained parallel execution, where each edge represents a task processed by a thread. Because threads in a warp execute in lockstep, distributing tasks and gathering results incurs minimal synchronization overhead. Unlike traditional approaches storing all child matches, this finegrained execution extends 𝑀 on the fly (expanding at most 32, the warp size) and avoids materializing all child matches at once and significantly reduces the space required for the execution stack. Figure 6 illustrates the fine-grained parallel execution model of gMatch. We use blue to denote the subtree rooted at 𝑣 57 and red for the subtree rooted at 𝑣 59 . While a warp attempts to extend the partial match 𝑀 = {(𝑢 1, 𝑣 58 ), (𝑢 2, 𝑣 57 )}, the intersection 𝑁 (𝑣 57 ) ∩ 𝑁 (𝑣 58 ) will be computed. As |𝑁 (𝑣 57 )| is smaller than |𝑁 (𝑣 58 )|, we say 𝐿 (𝑢 ) = 𝑁 (𝑣 ) and check the validity of each vertex within 𝐶𝑀 3 57 𝐿 (𝑢 ). Instead of computing and materializing all feasible can𝐶𝑀 3 didates, the warp only check the validity of the first 32 vertices

Shared Memory

Thread Block

Execution Stack of a Warp Buffers … … 𝑣33 𝑣34 … 𝑣50 𝑣51 … 𝑣56 … 𝑣60 𝑣61 …

𝑢1

𝑣58

𝑢2

𝑣57 𝑣59

𝑢3 𝑢4

32

Warp #0

Warp #1

Partial Matches

Warp #2

…

𝑣58 𝑣57 𝑣33 … 𝑣58 𝑣57 𝑣50 𝑣58 𝑣59 𝑣51

𝑁(𝑣57 ) 32

…

… …

𝑁(𝑣59 )

𝑣58 𝑣59 𝑣56

… …

𝑁(𝑣57 )

…

𝐿 (𝑢 ). All valid candidates among them (i.e., 𝑣 , 𝑣 , ..., 𝑣 ) are in 𝐶𝑀 3 1 2 32 stored in the buffer of the next level. The unexplored candidates in 𝐿 (𝑢 ) will be checked during backtracking. 𝐶𝑀 3 In summary, this design aligns naturally with the GPU’s architecture and offers two key benefits: (1) it reduces space overhead, as each thread handles only one candidate and does not need to store the feasible set or generate multiple partial matches; and (2) it enables fine-grained task divisibility and flexible scheduling across threads, the fundamental execution units on GPUs.

Warp #N-1

32

… …

…

Figure 7: Warp-level batch exploration.

𝜙 𝑢

𝐶

𝑣

𝑢

←

arg min ∈

𝑋

,

=𝑁 𝑣

𝑣

𝑢

𝑣 Output

𝑢

𝑣

𝑢 𝑣

𝑣

𝑣

𝑣

𝑣 … 𝑣 𝑣 𝑣 𝑣 × × 𝑣 𝑣 …

𝑣 ×

𝑣

𝑣

𝑣

𝑣

… 𝑣

𝑣 𝑣 …

𝑣

𝑣

× 𝑣

𝑣 ×

×

Figure 5: The search tree of fine-grained parallel execution on the example graph in Figure 1.

Shared Memory

Thread Block

Execution Stack of a Warp 𝑣58 𝑣57 𝑣59

𝑢3

𝑣1 𝑣2 𝑣3 𝑣4

𝑢4

Warp #1

… … … …

𝑁(𝑣57 ) 𝑣32

𝑣58 𝑣59

Warp #2

𝑁(𝑣59 )

𝑣1 𝑣2 … 𝑣32 … 𝑣50 𝑣58 𝑣60 𝑣61 𝑣51 𝑣52 … 𝑣56 𝑣58 𝑣61 𝑣62

…

…

𝑢1 𝑢2

Warp #0 Partial Matches 𝑣58 𝑣57

Buffers

Warp #N-1

Explored when backtracking

32

Figure 6: Fine-grained parallel execution on the example graphs in Figure 1. Remark. Existing GPU-based methods, such as STMatch, T-DFS, and G2 Miner, parallelize set intersection at the warp level using SIMT execution. Specifically, each thread processes one element from a set and checks its presence in the other set via binary search. This constitutes "fine-grained" parallelism within a single intersection. However, this strategy follows a coarse-grained execution model: a set intersection is treated as an indivisible unit, and all resulting candidates must be materialized in an array. As a result, each warp is assigned to extend a single partial match and enumerate all its feasible candidates, which leads to the issues described in Section 2.3. In contrast, our fine-grained computation model decomposes the extension of a partial match into independent percandidate validation tasks and assigns each task to a single thread. This allows a warp to simultaneously process tasks from multiple partial matches, or to distribute the processing of a single partial match across multiple warps, while significantly reducing memory overhead and improving the warp utilization.

4.2

Warp-Level Batch Exploration

Real-world graphs typically follow a power-law degree distribution, where most vertices have low degree. Moreover, for labeled graphs

and complex queries, advanced filtering techniques [34, 45] further reduce the number of candidates. As a result, the local candidate set of a partial match is often small, and its task group cannot fully utilize all threads in a warp, leading to underutilized GPU resources. Experimental results in Section 2.3 confirm the severity of this issue. To address this, we propose warp-level batch exploration on top of the fine-grained execution model. This model enables task divisibility and flexible scheduling across threads, making it natural to process multiple task groups in a batch. Instead of assigning a single task group to each warp, we unify multiple task groups into a shared virtual task pool. Each task group 𝑇𝑀𝐺 (𝑢) is defined 𝐿 (𝑢). During by a partial match 𝑀 and its local candidate set 𝐶𝑀 exploration, a warp maintains up to 32 partial matches at each level, with each thread responsible for one candidate mapping. These task groups collectively form a unified task pool accessible to all threads in the warp. Threads retrieve tasks from the pool using two lightweight pointers that track the current partial match and candidate position. This design incurs negligible memory overhead while significantly improving warp utilization without sacrificing memory efficiency. Materializing partial matches at each level can cause significant memory pressure on shared memory, with space complexity up to 𝑂 (|𝑉 (𝑄)| 2 ). To address this, we construct partial matches on-thefly by traversing the execution stack from bottom to top during the computation. Since each mapping from a query vertex to a feasible candidate defines a partial match, we store, for each candidate at level 𝑖 + 1, a pointer to its corresponding candidate at level 𝑖. This compact representation enables efficient reconstruction of partial matches with minimal memory overhead. Figure 7 illustrates the concept of warp-level batch exploration. This occurs during the backtracking phase of what is shown in Figure 6. To complete the tasks at level 2, the warp simultaneously processes the remaining part of 𝑁 (𝑣 57 ) and the entire 𝑁 (𝑣 59 ), treating each as a distinct batch. Since the number of unexplored candidates is fewer than 32, they can be processed within a single warp without requiring further backtracking. The feasible candidates for these two paths are stored at level 3 of the stack, denoted in blue and red, respectively. To extend these partial matches, the warp iterates through their local candidate sets using a fixed width of 32. For example, when processing the first 𝑁 (𝑣 57 ), the warp verifies the first 32 candidates and reports two final matches. The benefit of batch processing is more pronounced when processing 𝑁 (𝑣 59 ). To minimize thread idling, the warp batches candidates across multiple instances of 𝑁 (𝑣 59 ) together. This design inherently mitigates thread idling, as processing 𝑁 (𝑣 59 ) in a coarse-grained

manner would lead to severe underutilization. Importantly, such batch processing is not supported by prior coarse-grained models, where each partial match is treated as an indivisible task. Processing multiple partial matches together in those models would require significant memory for storing all candidate sets and intermediate states, e.g., the loop unrolling technique in STMatch [38].

4.3

Lightweight Load Balancing

The search space contains an exponential number of partial matches, each forming a task group of fine-grained tasks. This leads to an embarrassingly parallel workload. In our fine-grained execution model, intra-warp load balancing is naturally achieved through a unified task pool shared by threads within a warp. We therefore focus on load balancing across warps within a thread block and across thread blocks. Warps in the same thread block can coordinate via shared memory, while communication across blocks requires global memory. To fully utilize GPU resources, we design a lightweight two-phase load balancing mechanism: (1) an initialization phase that prepopulates an initial task pool to avoid early-stage underutilization, and (2) a dynamic balancing phase that redistributes work at runtime when imbalances occur. Initialization Phase. Although the total search space is large, the number of partial matches at the root level (i.e., first-level nodes in the search tree) may be small, leading to insufficient initial parallelism. To mitigate this, we first perform a breadth-first traversal of the search tree and accumulate partial matches into an initial task pool until the number of partial matches reaches a threshold 𝜏. We then launch the fine-grained parallel execution, where each warp fetches a partial match from the pool using an atomic counter when it runs out of work. The parameter 𝜏 is critical because it regulates the need for dynamic load balancing. When the number of initial partial matches significantly exceeds the number of available warps (typically hundreds to thousands on modern GPUs), load balancing becomes naturally easy because task variability is automatically smoothed out across a large shared task pool, reducing the necessity for explicit load-balancing mechanisms. This principle makes it straightforward to choose 𝜏 empirically: a larger value leads to better workload balance but comes with higher memory consumption. In our experiments, we set 𝜏 to 106 , achieving a good balance between memory consumption and efficiency. The experiments on the choice of 𝜏 are presented in Appendix B. Dynamic Balancing Phase. Once the initial task pool is exhausted, we shift to dynamic load redistribution. We follow the same highlevel idea as STMatch [38], idle warps receive half of the work from busy warps by splitting the execution stack of the selected busy warp, but our implementation is simpler and better aligned with our two-level load balancing. Specifically, when a warp runs out of tasks, it marks itself as idle in shared memory and enqueues its warp ID into a concurrent queue maintained per thread block. Each block also maintains a status array to track warp activity. Busy warps periodically check the queue before each recursive call. If a busy warp dequeues an idle warp ID, it splits its current task stack in half and assigns the split portion to the idle warp, updating the status accordingly. Concurrently, idle warps continuously poll their status in the array.

Once an idle warp detects its status change to active, it resumes execution with the newly assigned tasks. This avoids the costly stack scanning in STMatch and reduces contention, even though the busy warp selected may not always have the largest remaining workload. To support inter-block load balancing, the same mechanism is replicated at the block level using global memory for the queue and status array. In summary, by exploiting the intrinsic parallel structure of subgraph matching, our design makes work stealing a supplementary rather than primary mechanism for achieving load balancing.

5

IMPLEMENTATION AND COST ANALYSIS

In this section, we first introduce the implementation and then have a discussion on the cost of gMatch.

5.1

Implementation Details

Algorithm 2 presents the details of gMatch. We begin by performing a BFS to generate partial matches and populate the initial task pool, up to a threshold 𝜏. Each warp executes independently using an execution stack 𝑆, which is a 2D array of size |𝑉 (𝑄)| × 32. Each element 𝑆 [𝑖] [ 𝑗] represents the 𝑗-th thread’s state at depth 𝑖 and contains four fields: (1) 𝑣: a candidate data vertex, (2) 𝐹 : a flag indicating whether mapping 𝜙 [𝑖] to 𝑣 is valid, (3) 𝑝𝑖𝑑: the index of the parent vertex in 𝑆 [𝑖 − 1] corresponding to 𝜙 [𝑖 − 1], (4) 𝐶: a pointer to the local candidate set for the current partial match. Upon retrieving a partial match 𝑀 from the task pool, we initialize the execution stack 𝑆 and begin the search (Lines 2–4). Let 𝑙 denote the current recursion depth, where 𝜙 [𝑙] is the query vertex to be extended. At Line 7, we generate the virtual task pool. As shown in Lines 17–24, each thread traverses the execution stack to reconstruct its corresponding partial match on the fly without materializing it and computes the local candidate set for 𝜙 [𝑙]. Each resulting task group forms part of the virtual task pool, enabling fine-grained scheduling at the warp level. Each thread in the warp then retrieves a task from the virtual task pool (Lines 9–10). To avoid synchronization overhead, a designated leading thread assigns one task to each thread in the warp (Lines 25–32). Since a warp has only 32 threads, this assignment incurs negligible cost. Each thread then processes its assigned task by checking whether the candidate data vertex can extend the partial match. This involves verifying the label constraint, ensuring the vertex has not been used, and checking the edge constraints (Lines 33–41). If the generated partial match reaches the length of the matching order, the result is output (Lines 12–13). Otherwise, we use a warp-level collective primitive to check whether any thread has a feasible candidate; if so, the search continues (Lines 15–16).

5.2

Cost Analysis

We now analyze the space and time complexity of our parallel execution model and highlight its advantages over existing methods. Space Cost. Since all methods take the same input and produce the same output, we focus on the additional memory used during execution, specifically, the space required for the execution stack and load balancing structures. As the warp is the basic scheduling unit on GPUs, we analyze space usage at the warp level. Let 𝑊 be the number of threads per warp, and 𝐵 the number of warps per

Algorithm 2: gMatch Input : Query graph 𝑄, data graph 𝐺, matching order 𝜙, and initial task pool size 𝜏. Output : All matches from 𝑄 to 𝐺. 1 M ← generate 𝜏 partial matches with BFS-Search; /* Inter-warp parallelism */ 2 parallel foreach 𝑀 ∈ M do 3 Initialize execution stack 𝑆 of a warp with 𝑀; 4 Search(𝜙, 𝑆, |𝑀 |); Procedure Search(𝜙, 𝑆, 𝑙) 𝑖 ← 0, 𝑗 ← 0, 𝑘 ← 0; /* Each lane processes an element. 7 GenerateTask(𝜙, 𝑆, 𝑙, 𝐿𝐴𝑁 𝐸_𝐼 𝐷); 8 while 𝑡𝑟𝑢𝑒 do 9 if 𝐿𝐴𝑁 𝐸_𝐼 𝐷 = 0 then 10 (𝑖, 𝑗, 𝑘 ) ←ScatterTask(𝑆, 𝑙, 𝑖, 𝑗, 𝑘); 5

6

11 12 13 14 15 16

*/

/* Each lane processes a task. */ 𝑆 [𝑙 ] [𝐿𝐴𝑁 𝐸_𝐼 𝐷 ].𝐹 ←Process(𝜙, 𝑆, 𝑙, 𝑘, 𝐿𝐴𝑁 𝐸_𝐼 𝐷); if 𝑙 = |𝜙 | − 1 then if 𝑆 [𝑙 ] [𝐿𝐴𝑁 𝐸_𝐼 𝐷 ].𝐹 = 𝑡𝑟𝑢𝑒 then Output; else 𝐹 ← ballot_sync(0xFFFF FFFF, 𝑆 [𝑙 ] [𝐿𝐴𝑁 𝐸_𝐼 𝐷 ].𝐹 ); if 𝐹 ≠ 0 then Search(𝜙, 𝑆, 𝑙 + 1);

Procedure GenerateTask(𝜙, 𝑆, 𝑙, 𝑙𝑖𝑑) 𝑆 [𝑙 ] [𝑙𝑖𝑑 ].𝐶 ← ∅, 𝐶 ← ∅, 𝑝𝑖𝑑 ← 𝑙𝑖𝑑; 19 if 𝑆 [𝑙 − 1] [𝑙𝑖𝑑 ].𝐹 = 𝑓 𝑎𝑙𝑠𝑒 then return; 20 for 𝑖 ← 𝑙 − 1 to 0 do 21 𝑣 ← 𝑆 [𝑖 ] [𝑝𝑖𝑑 ].𝑣, 𝑝𝑖𝑑 ← 𝑆 [𝑖 ] [𝑝𝑖𝑑 ].𝑝𝑖𝑑;

17

18

𝜙

22 23 24

if 𝜙 [𝑖 ] ∈ 𝑁 + (𝜙 [𝑙 ] ) and (|𝑁 (𝑣) | ⩽ |𝐶 | or 𝐶 ≠ ∅) then 𝐶 ← 𝑁 (𝑣); 𝑆 [𝑙 ] [𝑙𝑖𝑑 ].𝐶 ← 𝐶;

Procedure ScatterTask(𝑆, 𝑙, 𝑖, 𝑗, 𝑘) while 𝑖 < 32 do 27 while 𝑗 < |𝑆 [𝑙 ] [𝑖 ].𝐶 | do 28 𝑆 [𝑙 ] [𝑘 ].𝑣 ← 𝑆 [𝑙 ] [𝑖 ].𝐶 [ 𝑗 ], 𝑆 [𝑙 ] [𝑘 ].𝑝𝑖𝑑 ← 𝑖; 29 𝑗 ← 𝑗 + 1, 𝑘 ← 𝑘 + 1; 30 if 𝑘 = 32 then return (𝑖, 𝑗, 𝑘 );

25

26

31 32

𝑖 ← 𝑖 + 1, 𝑗 ← 0; return (𝑖, 𝑗, 𝑘 );

Procedure Process(𝜙, 𝑆, 𝑙, 𝑘, 𝑙𝑖𝑑) if 𝑙𝑖𝑑 < 𝑘 then return 𝑓 𝑎𝑙𝑠𝑒; 35 𝑢 ← 𝜙 [𝑙 ], 𝑣 ← 𝑆 [𝑙 ] [𝑙𝑖𝑑 ].𝑣, 𝑝𝑖𝑑 ← 𝑆 [𝑙 ] [𝑙𝑖𝑑 ].𝑝𝑖𝑑; 36 if 𝐿 (𝑣) ≠ 𝐿 (𝑢 ) then return 𝑓 𝑎𝑙𝑠𝑒; 37 for 𝑖 ← 𝑙 − 1 to 0 do 38 𝑢 ′ ← 𝜙 [𝑖 ], 𝑣 ′ ← 𝑆 [𝑖 ] [𝑝𝑖𝑑 ].𝑣, 𝑝𝑖𝑑 ← 𝑆 [𝑖 ] [𝑝𝑖𝑑 ].𝑝𝑖𝑑;

33

34

if 𝑣 = 𝑣 ′ or (𝑢 ′ ∈ 𝑁 + (𝑢 ) and 𝑣 ∉ 𝑁 (𝑣 ′ )) then return 𝑓 𝑎𝑙𝑠𝑒; 𝜙

39 40 41

fits in shared memory. As a result, shared memory usage is minimal and does not limit warp occupancy. For load balancing, the local task queue within a thread block stores at most 𝐵 entries (one per warp), requiring 𝑂 (𝐵) space—also small enough for shared memory. The global task queue, used for inter-block coordination, is allocated in global memory. Its size is proportional to the number of blocks (typically hundreds), which is negligible compared to the available global memory (usually tens of gigabytes). Time Cost. The total time cost is determined by the size of the search space and the efficiency of processing each partial match. Given 𝑄, 𝐺, and 𝜙, the search space is fixed. We thus focus on analyzing the time to process a single partial match 𝑀. Let 𝑢 be the current query vertex to match. For each candidate 𝐿 (𝑢), feasibility checking involves two steps: (1) verifying that 𝑣 ∈ 𝐶𝑀 𝑣 is a neighbor of all data vertices mapped to 𝑢’s backward neigh𝜙 bors, i.e., for each 𝑢 ′ ∈ 𝑁 + (𝑢), we check whether 𝑣 ∈ 𝑁 (𝑀 [𝑢 ′ ])) via binary search, costing 𝑂 (log |𝑁 (𝑀 [𝑢 ′ ])|); (2) checking whether 𝑣 has already been used in 𝑀, which requires scanning the execution stack in 𝑂 (|𝑀 |) time. We fuse these checks into a single loop to improve efficiency. Thus, the total cost of processing 𝑀 is: Í Í 𝑂 ( 𝑣 ∈𝐶 𝐿 (𝑢 ) (|𝑀 | + 𝑢 ′ ∈𝑁 𝜙 (𝑢 ) log |𝑁 (𝑀 [𝑢 ′ ])|). 𝑀 + This cost is comparable to that of existing methods, but our approach achieves significantly better hardware utilization through fine-grained parallelism and warp-level batch exploration. The overhead of the scatter operation is negligible, as a warp processes at most 32 candidates. Advantages over Existing Methods. Compared with existing methods, our approach does not reduce the size of the search space or the time complexity of processing a partial match. Instead, it focuses on maximizing hardware utilization. First, our method reduces the space complexity of the execution stack per warp from 𝑂 (|𝑉 (𝑄)| ×𝑑 max ) to 𝑂 (|𝑉 (𝑄)|). This removes memory-related constraints on parallel execution, ensuring that the number of concurrent warps is not limited by stack size while also freeing substantial memory to handle larger graphs. Moreover, the stack can reside entirely in shared memory, improving access efficiency. Second, warp-level batch exploration enables full utilization of computing resources within a warp without incurring significant overhead, even when task groups are small. Finally, we employ a lightweight load balancing strategy that uses an initialization phase to prepopulate tasks, reducing the need for work stealing during execution.

6

EXPERIMENTS

In this section, we present the experimental settings and the results.

return 𝑡𝑟𝑢𝑒;

6.1 thread block. The execution stack depth equals |𝑉 (𝑄)|, and each thread maintains one candidate per level. Thus, the stack occupies 𝑂 (𝑊 × |𝑉 (𝑄)|) space per warp. Given that |𝑉 (𝑄)| is small (e.g., 12), the execution stack requires only a few kilobytes and comfortably

Experimental Setup

Methods Under Study. We compare gMatch with the latest GPUbased subgraph matching algorithms, including STMatch [38], EGSM [34], and T-DFS [43]. We also include two systems designed for mining small patterns: BEEP [23] and G2 Miner [8].

Their source code is publicly available on GitHub1,2,3,4,5 . To ensure fair comparison, we unify the matching order in STMatch and T-DFS by replacing their original strategies with the RI [6] ordering, consistent with our method. For EGSM, BEEP, and G2 Miner, we use their default ordering, which is specifically optimized for each framework. EGSM leverages a candidate graph for efficient candidate filtering and retrieval, which is particularly effective for workloads with large query graphs and medium-sized data graphs. We adopt the same candidate graph mechanism in gMatch when comparing against EGSM. However, we disable this feature when comparing against STMatch, T-DFS, BEEP, and G2 Miner, as they do not use candidate graphs and gain minimal benefit from candidate graphs for their target workloads. Our method’s compatibility with candidate graphs highlights its generality. Testbed. All competing algorithms are compiled with GCC 10.5 and NVCC 12.5. By default, we run our experiments on a Linux machine with an Intel Xeon Gold 6430 CPU and 256 GB main memory. The machine is equipped with an NVIDIA RTX 4090 GPU with 128 SMs and 24 GB global memory space. Datasets. We use 11 graph datasets in our experiments, as summarized in Table 3. These include eight real-world social networks (en, gh, gw, wt, lj, pk, fr, and ot), one co-authorship network (db), and two synthetic datasets (ld and rm). The dataset ld is from the LDBC social network benchmark [12]. We also use PaRMAT [25] to generate an RMAT graph rm because it can imitate the structure of real-world networks and allows us to control the graph scale to meet the memory constraints of a single GPU. The remaining datasets are from SNAP [28]. We choose the datasets en, gh, gw, wt, and db for comparison with EGSM, as they were also used in its experiments. These datasets are obtained directly from the EGSM authors’ release. Additionally, we include the datasets pk, fr, ot, and lj for comparison with STMatch and T-DFS, following the same datasets used in their original experiments. We choose the dataset ld as it is a popular benchmark for graph query and processing. We categorize the datasets based on their edge count: datasets with fewer edges belong to the first category, and those with more to the second. We use large query graphs on the first category and small query graphs on the second. All data graphs in the first category are labeled, where each vertex is assigned with a label uniformly at random. All datasets in the second category are unlabeled (i.e., |Σ| = 1).

P1

P5

P6

P4

P3

P2

P7

P8

P9

Figure 8: Query graphs. Query Sets. We evaluate the performance of subgraph matching algorithms across different query graph sizes. For small query graphs 1 https://github.com/HPC-Research-Lab/STMatch (commit: f98462a) 2 https://github.com/RapidsAtHKUST/EGSM (commit: de854d5) 3 https://github.com/lyuheng/tdfs (commit: 05ef1b1) 4 https://github.com/Nagi-Research-Group/BEEP (commit: 3eca170) 5 https://github.com/chenxuhao/GraphMiner (commit: 2a76e3f)

Table 3: The detailed statistics of datasets. Datasets

Abbr.

|𝑉 |

|𝐸|

|Σ|

𝑑 avg

𝑑 max

Enron GitHub Gowalla DBLP WikiTalk

en gh gw db wt

36,692 37,700 196,591 317,080 2,394,385

183,831 289,003 950,327 1,049,866 4,659,565

16 16 16 16 64

10.0 15.3 9.7 6.6 3.9

1,383 9,458 14,730 343 100,029

Pokec LiveJournal Orkut LDBC Friendster RMAT

pk lj ot ld fr rm

1,632,803 4,847,571 3,072,441 29,982,730 65,608,366 99,999,983

22,301,964 42,851,237 117,185,083 175,860,387 1,806,067,135 2,800,000,000

1 1 1 1 1 1

27.3 17.7 76.3 11.8 55.1 56.0

14,854 14,815 33,313 4,282,595 5,214 112,497

(|𝑉 (𝑄)| ≤ 5), we use the predefined query graphs shown in Figure 8. We consider 9 patterns which are commonly used in existing work [8, 19, 33, 43] for performance evaluation. All of them are unlabeled and will be used for testing on unlabeled data graphs. For larger query graphs, we generate 100 random queries per dataset using the same procedure as in EGSM and other prior work [5, 16, 33]. Specifically, to generate an 𝑛-vertex query graph, we start with a random seed vertex from the data graph and iteratively expand it by adding random neighboring vertices, along with all their connecting edges to the existing subgraph, until the target size of 𝑛 vertices is reached. Metrics. The total query time 𝑡 query for a query is divided into filtering time 𝑡 filtering , data transfer time (𝑡 transfer ), and search time (𝑡 search ). We define 𝑡 query as the time from when 𝐺 and 𝑄 are loaded into host memory until all results are found. 𝑡 filtering measures the time for candidate filtering and auxiliary structure construction, 𝑡 transfer captures the time to move 𝐺, 𝑄, and auxiliary data to GPU memory, and 𝑡 search accounts for the time to enumerate all matches. To ensure fair comparison, we report only 𝑡 search . This is because, for large queries, gMatch and EGSM share the same auxiliary structures, so any difference lies solely in 𝑡 search . For small queries, STMatch performs degree-based filtering on the CPU, which becomes a bottleneck on large datasets (e.g., ot), whereas gMatch and T-DFS implement the same filtering on the GPU, incurring negligible cost. Therefore, 𝑡 search provides a fair metric for comparison. To ensure all queries complete in a reasonable time, we set timeouts: 3 hours for small queries and 30 minutes for large queries. Queries exceeding the limit are marked as unsolved, and we report both the number of unsolved queries and the average search time. We also measure speedup using the following formula: 1 ∑︁ 𝑡 B (𝑄) Speedup(A over B) = , 𝑄 ∈ Q 𝑡 A (𝑄) |Q| where Q is the set of queries, and 𝑡 X (𝑄) is the search time of algorithm X on query 𝑄. We also evaluate GPU memory consumption. Since 𝐺 occupy the same amount of memory across all methods, we focus on measuring the space consumed by the execution stack.

6.2

Overall Comparison

Small Query. We evaluate gMatch against STMatch, T-DFS, BEEP and G2 Miner on the six large datasets. EGSM is excluded since it targets small data graphs and frequently runs out of memory

Table 4: Comparison of search time (millisecond) among T-DFS (TD), STMatch (ST), BEEP (BE), G2 Miner (G2), and gMatch (GM). OOT indicates timeouts, OOM indicates out-of-memory, and N/A indicates unsupported queries in the open-source code.

TD P1 P2 P3 P4 P5 P6 P7 P8 P9

G2

GM

TD

9,030 129 8,566 172 3,667 229 15,134 N/A 12,549 N/A 10,260 773 5,916 6,910 4,569 587 8,177 346

12 22 32 336 N/A N/A N/A N/A 55

23 86 323 1,626 4,390 7,935 3,295 1,646 288

112 2,083 5,665 12,932 926,807 3,865,506 842,156 520,329 106,059

G2

GM

TD OOM OOM OOM OOM OOM OOM OOM OOM OOM

TD

ST

ld BE

571 2,149 12,387 377,456 315,378 439,141 153,001 19,294 3,828

OOM OOM OOM OOM OOM OOM OOM OOM OOM

OOM 45 172 OOM 83 549 OOM OOM 5,844 N/A 205,991 152,834 N/A N/A 149,406 OOM N/A 167,778 OOM N/A 95,427 OOM N/A 10,191 OOM 114 910

en

gh

lj BE

ST

G2

GM

12,413 239 87 4,344 350 116 11,991 562 118 29,071 N/A 885 641,476 N/A N/A 509,852 16,875 N/A 259,700 23,703 N/A 153,237 12,218 N/A 20,799 3,995 4,347

31 422 2,080 4,127 334,373 637,667 305,968 221,493 25,114

ST

fr BE

G2

GM

OOM OOM OOM OOM OOM OOM OOM OOM OOM

OOM OOM OOM N/A N/A OOM OOM OOM OOM

1,224 4,354 4,090 12,132 OOT 53,556 OOT 2,806,445 N/A 1,068,031 N/A 1,519,540 N/A 731,884 N/A 235,798 4,643 33,776 gw

TD

ot BE

ST

G2

GM

463 12,486 661 92 5,601 9,587 840 190 19,244 43,329 1,380 299 OOT 285,055 N/A 9,730 1,234,177 4,153,864 N/A N/A 2,568,786 1,449,627 9,600 N/A 976,353 848,283 134,389 N/A 219,505 267,167 6,268 N/A 46,943 36,807 2,832 825

122 1,077 7,661 47,113 526,114 807,067 322,477 104,190 7,962

TD

ST

rm BE

G2

OOM OOM OOM OOM OOM OOM OOM OOM OOM

OOM OOM OOM OOM OOM OOM OOM OOM OOM

OOM OOM OOM N/A N/A OOM OOM OOM OOM

2,551 5,714 3,784 5,866 OOM 51,243 OOT 3,033,834 N/A 266,772 N/A 225,545 N/A 138,211 N/A 42,032 3,801 6,108

db

GM

wt

Time of gMatch (ms)

P1 P2 P3 P4 P5 P6 P7 P8 P9

84 431 987 5,473 13,429 29,846 13,225 6,223 1,754

pk BE

ST

Time of EGSM (ms)

Figure 9: Comparison of search time between EGSM and gMatch (millisecond). on these datasets. The result is shown in Table 4. For fairness, we disable auxiliary structures and BFS in gMatch, comparing only DFSbased search time with the baselines. gMatch is the only method that completes all queries across all datasets. In contrast, the baselines encounter OOM or OOT failures on large graphs (e.g., fr and rm) or graphs with high maximum degree (e.g., ld), demonstrating the superior scalability and generality of our approach. STMatch frequently runs out of memory because it allocates a buffer of size 𝑑 max at each stack level. In particular, it fails on ld, where 𝑑 max exceeds 4 million, and on fr and rm, where the remaining memory after stack allocation is insufficient. T-DFS is consistently slower than gMatch due to the overhead of its pagetable-based stack management. It also fails on fr and rm because representing data edges as arrays of edge pairs nearly doubles the memory footprint of the data graph, exceeding available memory. G2 Miner is a graph pattern mining system that supports different tasks such as motif discovery and frequent subgraph mining, with subgraph matching as a core operation. Given a query, it employs the AutoMine compiler [31] to optimize execution. However, the open-source release does not include the compiler and supports only a limited set of query patterns. As a result, our evaluation is

restricted to a subset of the patterns in Figure 8 (P1–P4 and P9). G2 Miner performs well on small graphs (e.g., pk, lj, and ot), but encounters OOM and OOT failures on large graphs. BEEP supports patterns with a single “central” vertex connected to all others (e.g., P1–P3 and P6–P9 in Figure 8) because it accelerates set intersection by constructing an adjacency matrix of size 2 ) for the neighborhood of the matched central vertex. This 𝑂 (𝑑 max matrix-based optimization enables good performance on P6–P9. However, it does not scale to large graphs and fails on ld, fr, and rm. Although the source code provides an option to disable this optimization, doing so leads to incorrect results. In summary, gMatch demonstrates superior scalability and broader applicability than existing methods. Large Query. To evaluate performance on larger query graphs, we conduct our experiments on en, gh, gw, db, and wt, as they are used in EGSM’s experiments. We compare only against EGSM, since the open-source implementations of the other methods do not support arbitrary queries with tens of vertices. For each dataset, we generate a query set consisting of 100 random 12-vertex query graphs using the method described in Section 6.1. Figure 9 demonstrates the result, where each red dot denotes a test case. Red dot below

the diagonal represents one test case where gMatch is faster than EGSM. Experiments result shows that gMatch outperforms EGSM on all datasets, achieving 36.58× average speedup on gw and 20.23× average speedup on db. This speedup mainly comes from improved warp utilization enabled by the fine-grained parallel execution and warp-level batch exploration. Memory Consumption. We compare gMatch’s execution-stack memory consumption with that of STMatch and T-DFS. Figure 10 shows the maximum memory consumption of three algorithms on each dataset, where each OOM denotes an out-of-memory error. We can see that STMatch incurs severe memory overhead. While T-DFS effectively reduces memory cost using a page table, it still consumes more memory than gMatch. gMatch minimizes the stack size to fit within GPU shared memory, requiring no additional global memory for the execution stack.

Figure 10: Memory consump- Figure 11: Speedup comparition of the execution stacks son vs. naive coarse-grained parallel execution. during DFS. Table 6: Comparison of idle rates between gMatch and the baseline STMatch under its different loop unrolling configurations (naive, unroll-2, unroll-4, unroll-8).

Table 5: Search time (millisecond) for all patterns across datasets. Each value represents the average search time over the lj, ot, pk, and fr datasets.

P1 P2 P3 P4 P5 P6 P7 P8 P9

Naive

Unroll-4

Unroll-8

Fine (naive)

Fine (batch)

1,355 4,334 20,284 935,433 617,268 924,871 418,832 186,048 21,460

1,277 4,098 19,022 877,175 578,856 866,658 392,459 174,420 20,106

1,221 3,917 18,443 849,963 560,772 839,577 380,196 168,945 19,487

1,249 4,058 18,855 869,191 573,618 858,787 388,898 174,327 19,922

1,132 3,429 15,905 714,827 483,227 743,052 340,906 140,781 16,785

Naive Unroll-2 Unroll-4 Unroll-8 gMatch

Evaluation of Individual Techniques

We evaluate the effectiveness of fine-grained parallel execution and warp-level batch exploration. We compare the performance of all query graphs from Figure 8. There are five methods under study: (1) Naive coarse-grained model; (2) Coarse-grained model with loop unrolling size of 4; (3) Coarse-grained model with loop unrolling size of 8; (4) Fine-grained model without batch exploration; and (5) Fine-grained model with batch exploration. We evaluate the performance of the five methods across six datasets: lj, ot, ld, pk, fr, and rm. The average search time for each query is reported in Table 5. Note that ld and rm are excluded from the average time calculation because the Naive, Unroll-4, and Unroll-8 methods fail on these datasets. The experimental results show that the fine-grained model with batch exploration consistently achieves the best performance. Only our fine-grained methods successfully handle the ld and rm datasets, which have high maximum degrees. All coarse-grained approaches fail on these challenging datasets due to out-of-memory errors. To evaluate the effectiveness of fine-grained parallel execution on large queries over medium-sized data graphs, we measure the average speedup on random query sets of 12 vertices (Figure 11). The baseline is the coarse-grained parallel approach without loop unrolling. Our batch exploration strategy achieves significant speedup

en

gw

gh

wt

70.74% 50.22% 30.16% 19.32% 4.14%

45.14% 31.38% 20.75% 10.84% 3.41%

40.34% 38.62% 22.94% 15.54% 3.53%

48.19% 59.47% 21.62% 11.41% 4.04%

58.80% 29.82% 12.49% 9.80% 1.51%

because candidate sets are small after filtering, causing substantial thread idling in coarse-grained models. To quantify this effect, we measure the idle rate (defined in Section 2.3) across the datasets, as shown in Table 6. Compared with existing methods, our approach reduces the idle rate to below 5%, demonstrating much higher efficiency.

6.4 6.3

db

Case Study

We conduct a case study on a production graph representing user payment relationships from our industry partner ByteDance, one of the largest social network companies in the world. In this graph, each vertex corresponds to a user, and each undirected edge denotes a payment transaction, reflecting mutual interaction. The graph contains 139,297,601 vertices and 153,853,965 edges. Figure 12a shows its degree distribution. Note that different from our previous experiments, we conduct all our case study on a V100 GPU with 32 GB device memory.

(a) Cumulative distribution of (b) Distribution of vertex risk lavertex degrees. bels.

Figure 12: Dataset of user payment relationships.

For small queries, we focus on identifying cycles of length 3 to 5, which represent circular transaction flows. These patterns are of particular interest as they often indicate suspicious or closed-loop behaviors. The detected cycles are used as input to downstream tasks such as graph neural network (GNN) analysis. The corresponding results are shown in Figure 13a. The graph also includes vertex labels indicating user risk levels, ranging from 0 to 14. Figure 12b illustrates the label distribution. For large-query evaluation, we analyze transaction patterns across different risk categories using queries with dozens of vertices. Due to privacy constraints, we cannot publish the actual queries and instead generate 30 synthetic queries with 12 vertices. We generate these queries by adapting the methodology from Section 6.1, with one modification: to ensure the inclusion of risky vertices, the seed vertex is always randomly selected from the set of vertices whose labels fall within the range of 10 to 14. Figure 13b presents the results of our case study. As shown in the figure, our method significantly outperforms the existing approaches.

(a) Small queries.

(b) Large queries.

Figure 13: Comparison of search times for the methods under study on the user payment network.

matching. While these works demonstrate high performance, they are not applicable to matching general patterns. CPU-based subgraph matching. Existing CPU-based subgraph matching algorithms primarily leverage pruning strategies and auxiliary data structures to enhance performance. The solutions of subgraph matching can be roughly divided into join-based and exploration-based, as indicated by [32]. Current exploration-based approaches are fundamentally rooted in Ullmann’s algorithms [37], incrementally matching query vertices to candidate vertices. Subsequent researches effectively reduce the search space by using optimized matching orders [6, 17], applying advanced filtering techniques [4, 26, 33, 41], avoid redundant computations [22, 29], and maintaining failure sets [3, 16]. While these exploration-based methods excel at handling complex queries on single machines through their complicated pruning strategies, join-based approaches offer better scalability for large-scale datasets by decomposing the query into sub-structures that can be processed in parallel, making them better suited for distributed environments. Graph pattern mining. Graph pattern mining (GPM) focuses on finding specific patterns within large data graphs. It includes several applications such as subgraph listing [8, 27], k-motif counting [8, 30], and k-frequent subgraph mining [8, 9]. Early work tackled these problems by creating specialized, custom solutions for each specific task. This strategy is effective but not general-purpose. To make it easier to run different mining tasks efficiently in parallel, generalpurpose graph mining systems were developed. These systems used an embedding-centric approach, starting with small pieces and gradually extending them step-by-step into larger subgraphs, checking at each stage if they meet the user’s criteria. Systems like [7, 35] used this model. However, similar to subgraph matching algorithms, these GPM systems also face the problem of tracking the huge number of intermediate partial results. To solve this problem, [8, 10, 31] adopt DFS strategies for enumeration.

8 7

RELATED WORK

GPU-based subgraph matching. Recent research on subgraph matching has increasingly leveraged GPUs due to their massive parallelism capabilities. Existing GPU-based subgraph matching solutions can be categorized into BFS-based, DFS-based, and hybrid approaches according to their search strategy. Early GPU subgraph matching research primarily adopted BFS-based methods due to their inherent workload balance advantages [36, 40, 45]. However, BFS-based approaches typically generate a large amount of intermediate results that cannot fit in the memory of GPUs. To reduce memory consumption, recent studies have revisited DFS-based strategies [11, 23, 38, 43] or DFS-BFS hybrid strategies [34]. Moreover, there are some researches focusing on processing subgraph matching on very large data graphs that cannot fit in GPU memory. They either divide the data graph into multiple partitions [14] or use a view-based strategy according to the matched data vertex of the source query vertex [20, 42, 44]. These algorithms preprocess the data graph by dividing the problem into multiple chunks and load one chunk into memory at a time. Additionally, there are GPU algorithms specifically designed for 𝑘-clique counting [1, 24] and triangle counting [2, 13, 18], which are special cases of subgraph

CONCLUSION

We presented gMatch, a fine-grained and hardware-efficient subgraph matching method for GPUs. Our work is driven by the key insight that the coarse-grained parallel execution model fundamentally conflicts with real-world power-law graphs, giving rise to two inherent issues. First, GPU scalability is constrained by memory: large per-task stacks severely limit the number of concurrent warps and underutilize GPU parallelism, making a lightweight stack design essential. Second, mapping highly irregular graph workloads onto fixed-size warps leads to significant thread underutilization. To address these challenges, our fine-grained design incorporates warp-level batch exploration and a lightweight loadbalancing mechanism, achieving high efficiency and scalability across a wide range of workloads. Experimental results show that gMatch delivers strong performance and scalability with low memory overhead, making it well suited for practical, large-scale graph applications.

ACKNOWLEDGMENTS The work is supported by the National Natural Science Foundation of China (NSFC) under Grant No. 62572303. Shixuan Sun and Minyi Guo are the corresponding authors.

REFERENCES [1] Mohammad Almasri, Izzat El Hajj, Rakesh Nagi, Jinjun Xiong, and Wen-mei Hwu. 2022. Parallel K-clique counting on GPUs. In Proceedings of the 36th ACM International Conference on Supercomputing (Virtual Event) (ICS ’22). Association for Computing Machinery, New York, NY, USA, Article 21, 14 pages. [2] Mohammad Almasri, Neo Vasudeva, Rakesh Nagi, Jinjun Xiong, and Wen-Mei Hwu. 2021. HyKernel: A Hybrid Selection of One/Two-Phase Kernels for Triangle Counting on GPUs. In 2021 IEEE High Performance Extreme Computing Conference (HPEC). 1–7. [3] Junya Arai, Yasuhiro Fujiwara, and Makoto Onizuka. 2023. GuP: Fast Subgraph Matching by Guard-based Pruning. Proc. ACM Manag. Data 1, 2, Article 167 (June 2023), 26 pages. [4] Bibek Bhattarai, Hang Liu, and H. Howie Huang. 2019. CECI: Compact Embedding Cluster Index for Scalable Subgraph Matching. In Proceedings of the 2019 International Conference on Management of Data (Amsterdam, Netherlands) (SIGMOD ’19). Association for Computing Machinery, New York, NY, USA, 1447–1462. [5] Fei Bi, Lijun Chang, Xuemin Lin, Lu Qin, and Wenjie Zhang. 2016. Efficient Subgraph Matching by Postponing Cartesian Products. In Proceedings of the 2016 International Conference on Management of Data, SIGMOD Conference 2016, San Francisco, CA, USA, June 26 - July 01, 2016, Fatma Özcan, Georgia Koutrika, and Sam Madden (Eds.). ACM, 1199–1214. [6] Vincenzo Bonnici, Rosalba Giugno, Alfredo Pulvirenti, Dennis Shasha, and Alfredo Ferro. 2013. A subgraph isomorphism algorithm and its application to biochemical data. BMC bioinformatics 14 (2013), 1–13. [7] Hongzhi Chen, Miao Liu, Yunjian Zhao, Xiao Yan, Da Yan, and James Cheng. 2018. G-Miner: an efficient task-oriented graph mining system. In Proceedings of the Thirteenth EuroSys Conference (Porto, Portugal) (EuroSys ’18). Association for Computing Machinery, New York, NY, USA, Article 32, 12 pages. [8] Xuhao Chen and Arvind. 2022. Efficient and Scalable Graph Pattern Mining on GPUs. In 16th USENIX Symposium on Operating Systems Design and Implementation (OSDI 22). USENIX Association, Carlsbad, CA, 857–877. [9] Xuhao Chen, Roshan Dathathri, Gurbinder Gill, and Keshav Pingali. 2020. Pangolin: an efficient and flexible graph mining system on CPU and GPU. Proc. VLDB Endow. 13, 8 (April 2020), 1190–1205. [10] Vinicius Dias, Carlos H. C. Teixeira, Dorgival Guedes, Wagner Meira, and Srinivasan Parthasarathy. 2019. Fractal: A General-Purpose Graph Pattern Mining System. In Proceedings of the 2019 International Conference on Management of Data (Amsterdam, Netherlands) (SIGMOD ’19). Association for Computing Machinery, New York, NY, USA, 1357–1374. [11] Vibhor Dodeja, Mohammad Almasri, Rakesh Nagi, Jinjun Xiong, and Wen-mei Hwu. 2022. PARSEC: PARallel Subgraph Enumeration in CUDA. In 2022 IEEE International Parallel and Distributed Processing Symposium (IPDPS). 168–178. [12] Orri Erling, Alex Averbuch, Josep-Lluis Larriba-Pey, Hassan Chafi, Andrey Gubichev, Arnau Prat-Pérez, Minh-Duc Pham, and Peter A. Boncz. 2015. The LDBC Social Network Benchmark: Interactive Workload. In SIGMOD. 619–630. [13] Oded Green, Pavan Yalamanchili, and Lluís-Miquel Munguía. 2014. Fast Triangle Counting on the GPU. In 2014 4th Workshop on Irregular Applications: Architectures and Algorithms (IA3 ). 1–8. [14] Wentian Guo, Yuchen Li, Mo Sha, Bingsheng He, Xiaokui Xiao, and Kian-Lee Tan. 2020. GPU-Accelerated Subgraph Enumeration on Partitioned Graphs. In Proceedings of the 2020 ACM SIGMOD International Conference on Management of Data (Portland, OR, USA) (SIGMOD ’20). Association for Computing Machinery, New York, NY, USA, 1067–1082. [15] Pankaj Gupta, Venu Satuluri, Ajeet Grewal, Siva Gurumurthy, Volodymyr Zhabiuk, Quannan Li, and Jimmy Lin. 2014. Real-time twitter recommendation: Online motif detection in large dynamic graphs. Proceedings of the VLDB Endowment 7, 13 (2014), 1379–1380. [16] Myoungji Han, Hyunjoon Kim, Geonmo Gu, Kunsoo Park, and Wook-Shin Han. 2019. Efficient Subgraph Matching: Harmonizing Dynamic Programming, Adaptive Matching Order, and Failing Set Together. In SIGMOD Conference. 1429–1446. [17] Huahai He and Ambuj K. Singh. 2008. Graphs-at-a-time: query language and access methods for graph databases. In Proceedings of the 2008 ACM SIGMOD International Conference on Management of Data (Vancouver, Canada) (SIGMOD ’08). Association for Computing Machinery, New York, NY, USA, 405–418. [18] Yang Hu, Hang Liu, and H. Howie Huang. 2018. High-Performance Triangle Counting on GPUs. In 2018 IEEE High Performance extreme Computing Conference (HPEC). 1–5. [19] Kasra Jamshidi, Rakesh Mahadasa, and Keval Vora. 2020. Peregrine: a patternaware graph mining system. In Proceedings of the Fifteenth European Conference on Computer Systems (Heraklion, Greece) (EuroSys ’20). Association for Computing Machinery, New York, NY, USA, Article 13, 16 pages. [20] Guanxian Jiang, Qihui Zhou, Tatiana Jin, Boyang Li, Yunjian Zhao, Yichao Li, and James Cheng. 2022. VSGM: view-based GPU-accelerated subgraph matching on large graphs. In Proceedings of the International Conference on High Performance Computing, Networking, Storage and Analysis (Dallas, Texas) (SC ’22). IEEE Press, Article 52, 15 pages.

[21] Jiawei Jiang, Yusong Hu, Xiaosen Li, Wen Ouyang, Zhitao Wang, Fangcheng Fu, and Bin Cui. 2022. Analyzing Online Transaction Networks with Network Motifs. In Proceedings of the 28th ACM SIGKDD Conference on Knowledge Discovery and Data Mining (Washington DC, USA) (KDD ’22). Association for Computing Machinery, New York, NY, USA, 3098–3106. [22] Tatiana Jin, Boyang Li, Yichao Li, Qihui Zhou, Qianli Ma, Yunjian Zhao, Hongzhi Chen, and James Cheng. 2023. Circinus: Fast Redundancy-Reduced Subgraph Matching. Proc. ACM Manag. Data 1, 1, Article 12 (May 2023), 26 pages. [23] Samiran Kawtikwar, Mohammad Almasri, Wen-Mei Hwu, Rakesh Nagi, and Jinjun Xiong. 2023. Beep: Balanced efficient subgraph enumeration in parallel. In Proceedings of the 52nd International Conference on Parallel Processing. 142–152. [24] Vinayak Kesarwani, Shivangi Gaur, Sudeep Ranjan Sahoo, Kishan Tamboli, and Vishwesh Jatala. 2024. A Partitioning Scheme for Large Scale Clique Counting on Single GPU. In 2024 IEEE 31st International Conference on High Performance Computing, Data and Analytics Workshop (HiPCW). 167–168. [25] Farzad Khorasani, Rajiv Gupta, and Laxmi N. Bhuyan. 2015. Scalable SIMDEfficient Graph Processing on GPUs. In Proceedings of the 24th International Conference on Parallel Architectures and Compilation Techniques (PACT ’15). 39– 50. [26] Hyunjoon Kim, Yunyoung Choi, Kunsoo Park, Xuemin Lin, Seok-Hee Hong, and Wook-Shin Han. 2021. Versatile Equivalences: Speeding up Subgraph Query Processing and Subgraph Matching. In Proceedings of the 2021 International Conference on Management of Data (Virtual Event, China) (SIGMOD ’21). Association for Computing Machinery, New York, NY, USA, 925–937. [27] Hyeonji Kim, Juneyoung Lee, Sourav S. Bhowmick, Wook-Shin Han, JeongHoon Lee, Seongyun Ko, and Moath H.A. Jarrah. 2016. DUALSIM: Parallel Subgraph Enumeration in a Massive Graph on a Single Machine. In Proceedings of the 2016 International Conference on Management of Data (San Francisco, California, USA) (SIGMOD ’16). Association for Computing Machinery, New York, NY, USA, 1231–1245. [28] Jure Leskovec and Andrej Krevl. 2014. SNAP Datasets: Stanford Large Network Dataset Collection. http://snap.stanford.edu/data. Accessed: April 11, 2026. S Subgraph Matching [29] Yujie Lu, Zhijie Zhang, and Weiguo Zheng. 2025. B○X: with Batch Backtracking Search. Proc. ACM Manag. Data 3, 1, Article 15 (Feb. 2025), 27 pages. [30] Chenhao Ma, Reynold Cheng, Laks V. S. Lakshmanan, Tobias Grubenmann, Yixiang Fang, and Xiaodong Li. 2019. LINC: a motif counting algorithm for uncertain graphs. Proc. VLDB Endow. 13, 2 (Oct. 2019), 155–168. [31] Daniel Mawhirter and Bo Wu. 2019. AutoMine: harmonizing high-level abstraction and high performance for graph mining. In Proceedings of the 27th ACM Symposium on Operating Systems Principles (Huntsville, Ontario, Canada) (SOSP ’19). Association for Computing Machinery, New York, NY, USA, 509–523. [32] Shixuan Sun and Qiong Luo. 2020. In-Memory Subgraph Matching: An Indepth Study. In Proceedings of the 2020 ACM SIGMOD International Conference on Management of Data (Portland, OR, USA) (SIGMOD ’20). Association for Computing Machinery, New York, NY, USA, 1083–1098. [33] Shixuan Sun, Xibo Sun, Yulin Che, Qiong Luo, and Bingsheng He. 2020. RapidMatch: a holistic approach to subgraph query processing. Proc. VLDB Endow. 14, 2 (Oct. 2020), 176–188. [34] Xibo Sun and Qiong Luo. 2023. Efficient GPU-Accelerated Subgraph Matching. Proc. ACM Manag. Data 1, 2, Article 181 (June 2023), 26 pages. [35] Carlos H. C. Teixeira, Alexandre J. Fonseca, Marco Serafini, Georgos Siganos, Mohammed J. Zaki, and Ashraf Aboulnaga. 2015. Arabesque: a system for distributed graph mining. In Proceedings of the 25th Symposium on Operating Systems Principles (Monterey, California) (SOSP ’15). Association for Computing Machinery, New York, NY, USA, 425–440. [36] Ha-Nguyen Tran, Jung-jae Kim, and Bingsheng He. 2015. Fast Subgraph Matching on Large Graphs using Graphics Processors. In Database Systems for Advanced Applications, Matthias Renz, Cyrus Shahabi, Xiaofang Zhou, and Muhammad Aamir Cheema (Eds.). Springer International Publishing, Cham, 299–315. [37] J. R. Ullmann. 1976. An Algorithm for Subgraph Isomorphism. J. ACM 23, 1 (Jan. 1976), 31–42. [38] Yihua Wei and Peng Jiang. 2022. STMatch: accelerating graph pattern matching on GPU with stack-based loop optimizations. In Proceedings of the International Conference on High Performance Computing, Networking, Storage and Analysis (Dallas, Texas) (SC ’22). IEEE Press, Article 53, 13 pages. [39] Peter Willett. 1999. Matching of Chemical and Biological Structures Using Subgraph and Maximal Common Subgraph Isomorphism Algorithms. Springer New York, New York, NY, 11–38. [40] Lizhi Xiang, Arif Khan, Edoardo Serra, Mahantesh Halappanavar, and Aravind Sukumaran-Rajam. 2021. cuTS: scaling subgraph isomorphism on distributed multi-GPU systems using trie based data structure. In Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis (St. Louis, Missouri) (SC ’21). Association for Computing Machinery, New York, NY, USA, Article 69, 14 pages. [41] Rongjian Yang, Zhijie Zhang, Weiguo Zheng, and Jeffrey Xu Yu. 2023. Fast Continuous Subgraph Matching over Streaming Graphs via Backtracking Reduction. Proc. ACM Manag. Data 1, 1, Article 15 (May 2023), 26 pages.

[42] Lyuheng Yuan, Akhlaque Ahmad, Da Yan, Jiao Han, Saugat Adhikari, Xiaodong Yu, and Yang Zhou. 2024. G2-AIMD: A Memory-Efficient Subgraph-Centric Framework for Efficient Subgraph Finding on GPUs. In 2024 IEEE 40th International Conference on Data Engineering (ICDE). 3164–3177. [43] Lyuheng Yuan, Da Yan, Jiao Han, Akhlaque Ahmad, Yang Zhou, and Zhe Jiang. 2024. Faster Depth-First Subgraph Matching on GPUs. In 2024 IEEE 40th International Conference on Data Engineering (ICDE). 3151–3163. [44] Li Zeng, Lei Zou, and M. Tamer Özsu. 2023. SGSI – A Scalable GPU-Friendly Subgraph Isomorphism Algorithm. IEEE Transactions on Knowledge and Data Engineering 35, 11 (2023), 11899–11916.

[45] Li Zeng, Lei Zou, M. Tamer Özsu, Lin Hu, and Fan Zhang. 2020. GSI: GPUfriendly Subgraph Isomorphism. In 2020 IEEE 36th International Conference on Data Engineering (ICDE). 1249–1260. [46] Shijie Zhang, Shirong Li, and Jiong Yang. 2009. GADDI: distance index based subgraph matching in biological networks. In Proceedings of the 12th International Conference on Extending Database Technology: Advances in Database Technology (Saint Petersburg, Russia) (EDBT ’09). Association for Computing Machinery, New York, NY, USA, 192–203. [47] Zhijie Zhang, Yujie Lu, Weiguo Zheng, and Xuemin Lin. 2024. A Comprehensive Survey and Experimental Study of Subgraph Matching: Trends, Unbiasedness, and Interaction. Proc. ACM Manag. Data 2, 1, Article 60 (March 2024), 29 pages.

EGSM gMatch

105 104 103 102 101

Dense Sparse db

Dense Sparse gw

107 106 105 104 103 102 101

EGSM gMatch

8

12 db

16

8

12 gw

16

100

EGSM gMatch

80 60 40 20 0

8

12 db

16

8

12 gw

16

Average Search Time (ms)

(a) Search time with varying 𝑄 (b) Search time with varying |𝑉 (𝑄 ) |. density. 107 106 105 104 103 102 101

EGSM gMatch

12

16 db

20

12

16 gw

20

100 #Solved Queries

EGSM gMatch

80 60 40 20 0

12

16 db

20

12

16 gw

20

Average Search Time (ms)

(c) Number of solved queries with (d) Search time with varying varying |𝑉 (𝑄 ) |. |Σ(𝐺 ) |.

𝑖=1 𝑖

label (an integer within [0,|Σ| −1]) and 𝛼 controls skewness. We test 𝛼 = {0.6, 0.8, 1.0}, generating progressively skewed distributions. Specifically, a larger 𝛼 denotes more skewed distribution. To ensure that all queries can be solved within the time limit, we execute 100 random query graphs with 10 vertices for each configuration. The result is shown in Figure 14f. Analysis of Data Graph Sizes. To evaluate how gMatch and competing methods perform as data graph size increases, we use the LDBC benchmark to generate unlabeled, undirected data graphs of varying sizes, using LDBC’s standard scale factors of 3, 10, 30, and 100. These parameters directly determine the graph size and were chosen to ensure that all graphs fit within the memory of a single GPU, as all evaluated methods operate in an in-memory setting. The key statistics of the generated graphs are summarized in Table 7.

106

Average Search Time (ms)

Given that the large 𝑑 max in the LDBC dataset causes STMatch to run out of memory, we exclude it and only evaluate gMatch and T-DFS in this experiment. Similar to Section 6.2, we use the same matching order and symmetry breaking strategy in gMatch and T-DFS for fairness of comparison. As shown in Figure 15, gMatch demonstrates robust scalability across data graphs of varying sizes. In contrast, T-DFS is slower than gMatch in all cases and, moreover, failed to complete execution for all patterns at scale factor 100 and for pattern P4 at scale factors 3 and 30. This failure is due to excessive memory allocation by its memory management system. These results collectively confirm the efficiency and scalability of gMatch.

Average Search Time (ms)

SCALABILITY EVALUATIONS

#Solved Queries

A

To evaluate the scalability of our algorithm, we conduct four experiments on two representative datasets (db and gw) and one experiment using the LDBC benchmark. We consider five aspects: query graph density, query graph size, label cardinality, skewness of data vertex label, and data graph size. For analysis of data graph sizes, we conduct this experiment exclusively with T-DFS, as EGSM fails to process large data graph and STMatch does not support large maximum degree. We conduct the rest of scalability evaluations exclusively with EGSM, since the released code of STMatch and T-DFS only support small query graphs with no more than 7 vertices. Analysis of Query Graph Density. We first examine how query graph density impacts performance. Using 100 query graphs with 12 vertices per dataset, we categorize query graphs based on average degree (𝑑 avg ): sparse (𝑑 avg < 3) and dense (𝑑 avg ≥ 3). As shown in Figure 14a, the query time is generally longer with smaller query graph density. This is because sparse query graphs usually lead to larger intermediate results and are more difficult to solve. Analysis of Query Graph Size. To evaluate the impact of query graph size, we generate 100 random query graphs per dataset at each of three sizes: 8, 12 and 16 vertices. The average query time of solved queries is shown in Figure 14b, and the number of solved queries is shown in Figure 14c. Note that all queries with size 8 or 12 can be solved within the time limit. Typically, a larger query graph is more difficult to solve than smaller ones. Analysis of Label Cardinality. To examine how the number of distinct labels |Σ(𝐺)| affects performance, we use the configurations with |Σ(𝐺)| = 12, 16, and 20. For each configuration, we generate 100 random 12-vertex queries. The average query time of solved queries is shown in Figure 14d, and Figure 14e shows the number of solved queries. Note that all queries under |Σ(𝐺)| = 16, 20 are solved. Generally speaking, a larger label cardinality leads to simpler subgraph matching, because this leads to smaller candidates set and thus increases the selectivity. Analysis of Label Distribution Skewness. We evaluate how skewed label distributions impact performance by synthetically modifying labels of db and gw using Zipf’s law. For each dataset, we first fix the label set size |Σ(𝐺)| = 16, then assign vertex labels −𝛼 according to this distribution: 𝑃 (label = ℓ) = Í(ℓ+1) |Σ| −𝛼 , where ℓ is the

107 106 105 104 103 102 101

EGSM gMatch

0.6 0.8 1.0 db

0.6 0.8 1.0 gw

(e) Number of solved queries with (f) Search time with varying data varying |Σ(𝑄 ) |. graph skewness.

Figure 14: Scalability evaluations.

|𝑉 |

|𝐸|

𝑑 avg

𝑑 max

3 10 30 100

9,281,922 29,982,730 88,789,833 282,386,021

52,651,300 175,860,387 540,506,176 1,773,425,640

11.34 11.73 12.18 12.56

1,346,287 4,282,595 12,684,688 40,767,884

Related documents

• Property Testing for Recursive Query Languages#639749
• From Data Querying to Data Investigations: Rethinking Natural Language Interfaces for Databases#639750
• Aker: Density-Aware Approximate Caching for Vector Search (Extended Version)#639751
• Convolution Sum-Product Queries#639752
• Tree Databases#639753