HDReason: Algorithm-Hardware Codesign for Hyperdimensional Knowledge Graph Reasoning

Figure 2(a) illustrates the basic concepts of KG and provides an example of the KGC task. A KG can be formally represented as a directed graph $G=\{(v,\,r,\,u)\>|\>u,v\in\xi,\>r\in R\}$, where $\xi$ is the set of entities (vertices) and $R$ is the set of relations (vashishth2019composition, ). Each directed edge in a KG represents a factual statement and is defined as a fact triple $l=(v,\,r,\,u)$. As shown in Figure 2, (Anne Hathaway, born_in, New York City) is a fact triple, which presents the statement Anne Hathaway was born in the New York City. Notice that relations are directed edges to distinguish the subject and object in factual statements.

In various KG-related reasoning applications, we focus on the KGC task due to its elevated computational demands and training challenges (chen2020review, ). The KGC task is also referred to as KG link prediction. As a low-level logic reasoning (gao2021quatde, ), KGC aims to infer nonexistent relations from available knowledge. An example task is to find a new triple like (Anne Hathaway, born_in, ?) in the context of Figure 2(a).

Prior works generally handle two types of KGC: single direction reasoning and double direction reasoning. The former is limited to finding a possible connection from the source vertex to the destination vertex. And the latter requires finding both positive connection (v,r,?) and negative connection (?,r,u).

The score function measures the distance between two nodes given its relation type (zheng2020dgl, ) in a KG link prediction task. For a subject $\vec{e^{v}_{i}}$ and the relation $\vec{e^{r}_{k}}$, if vertex $j$ is connected with vertex $i$ via relation $r_{k}$, the mapping function f is defined as:

(3)

$$\vec{e^{v}_{j}}\approx f(\vec{e^{v}_{i}},\vec{e^{r}_{k}})$$

Figure 1(d) provides an example of TransE (bordes2013translating, ), a type of score function. TransE treats the relationship as a translation vector so that the embedded entities are connected by relation r with low error. The TransE score function is defined as:

(4)

$$score_{j}=N_{x}(\vec{e^{v}_{i}}+\vec{e^{r}_{k}}-\vec{e^{v}_{j}})$$

Here $N_{x}$ stands for the norm function. $\vec{e^{v}_{i}}$, $\vec{e^{r}_{k}}$, and $\vec{e^{v}_{j}}$ can be also treated as subject vector, relation vector, and object vector. From a semantic point of view, equations 3 and 4 stand for the subject being connected with the object via a specific relation. Given the query subject and relation, we use the TransE score function to calculate the score of every vertex in the graph. A larger score means that this vertex has a higher chance of being connected with the query subject vertex via the query relation.

In table 1, we summarize prior efforts that try to speed up GNN models. Historically, most GNN-based model accelerations leveraged GPUs built on the PyG framework (vashishth2019composition, ). While GPUs are versatile computing platforms, they suffer from static data paths and reduced energy efficiency (lacey2016deep, ). In contrast, FPGAs offer a balance between computational parallelism and energy consumption, making them suitable for enhancing GNN-based models. Nevertheless, it is challenging to apply existing FPGA-based GCN training accelerators directly to solve KGC tasks (zeng2020graphact, ; lin2022hp, ). One major challenge is that KG datasets incorporate semantic data for both vertices and edges. The cutting-edge FPGA-based accelerator, FlowGNN (sarkar2023flowgnn, ), does support these embeddings, but solely for GNN model inference. Furthermore, current FPGA accelerators typically focus on GNN propagation weight training, without supporting the crucial vertex and relation embedding training. Given the semantic richness of KG datasets, this training is central to the entire reasoning process. KGC training on GPU suffers from large-scale GPU memory usage for updating large-scale embeddings, extensive epochs for learning semantic information, and complexities in integrating accelerated score functions and graph models on FPGA platforms. A hardware-software approach is vital for an optimal FPGA framework tailored for KGC training.

The embedding vectors for all vertices comprise the embedding matrix $\textbf{e}^{\textbf{v}}$, with a dimension of $|V|\times d$. Similarly for the relations, we have an embedding matrix $\textbf{e}^{\textbf{r}}$ of size $|R|\times d$. To map all vertices and relations from the original embedding space into hyperspace, we have the matrix multiplication equation:

(5)

$$\mathbf{H^{v}}=tanh(\mathbf{e^{v}}\mathbf{H^{B}})$$

(6)

$$\mathbf{H^{r}}=tanh(\mathbf{e^{r}}\mathbf{H^{B}})$$

$\textbf{H}^{\textbf{v}}$ and $\textbf{H}^{\textbf{r}}$ are the encoded vertex and relation hypervector matrix with dimension $|V|\times D$ and $|R|\times D$ respectively.

After generating the vertex and relation hypervector matrix, we use an aggregation operation over the vertex and its neighbors to generate the corresponding memory hypervector as stated in equation 1. Given that each edge now also has a hyperspace mapping, it’s essential to first execute binding operations between each vertex hypervector and its related relation hypervectors, ensuring the vertex data is correctly linked with its relevant relation information:

(7)

$$\vec{M^{v}_{i}}=\Sigma_{(j,r)\in N(i)}\vec{H^{v}_{j}}\circ\vec{H^{r}_{r}}$$

, where $\circ$ represents element-wise Hadamard product. Equation 7 can also be rewritten in a matrix-to-matrix multiplication (M2MM) format:

(8)

$$\mathbf{M^{v}}=\Sigma_{r\in R}((\mathbf{A^{r}}\mathbf{H^{v}})\circ\mathbf{E^{r}})$$

$\mathbf{A^{r}}$ is the adjacency matrix over relation r, with a dimension of $|V|\times|V|$. If there exists an entity ($v_{i}$, $r$, $v_{j}$) $\in$ $\xi$, then $A^{r}_{ij}=1$, otherwise $A^{r}_{ij}=0$. Matrix $\mathbf{E^{r}}$ with dimension $|V|\times D$ is the concatenation of relation embedding hypervectors across all vertices:

(9)

$$\vec{E^{r}_{i}}=\vec{e_{r}}\quad i\in[0:|V|-1]$$

After obtaining vertex memory hypervectors, we use a score function such as the backend decoder to perform link prediction tasks. This encoder-decoder structure is widely adopted by previous GNN-based works (schlichtkrull2018modeling, ; shang2019end, ; vashishth2019composition, ). Here we choose to use TransE (bordes2013translating, ) as the score function. Assuming the input query is subject $\vec{e^{v}_{i}}$ and relation $\vec{e^{r}_{k}}$. To find the object vertex that is connected to subject vertex $v^{i}$ via relation $r^{k}$, we use the following vertex score function:

(10)

$$\vec{P}=Sigmoid(|M^{v}_{i}+H^{r}_{k}-\mathbf{M^{v}}|_{1}+bias)$$

The score vector P has a shape of $|V|\times 1$, with its element representing the possibility of a vertex connecting with vertex $v_{i}$ through relation $r$. In the equation above, the vectors $M^{v}_{i}$ and $H^{r}_{k}$ are first broadcasted to the same shape as the matrix $\mathbf{M^{v}}$. The L1 normalization is carried out only on one axis of the matrix, giving a vector form.

On top of previous HDC graph learning works (poduval2022graphd, ; nunes2022graphhd, ; kang2022relhd, ), we introduce the original space embedding to vertices and relations; we also add a score function to the HDC graph model. Now we can conduct backpropagation and train the original embedding model ($\textbf{e}^{\textbf{v}}$ and $\textbf{e}^{\textbf{r}}$):

(11)

$$\nabla_{\mathbf{e^{v}}}\mathbf{L}=\frac{\partial\mathbf{L}}{\partial\mathbf{P}}\frac{\partial\mathbf{P}}{\partial\mathbf{M^{v}}}\frac{\partial\mathbf{M^{v}}}{\partial\mathbf{H^{v}}}\frac{\partial\mathbf{H^{v}}}{\partial\mathbf{e}^{v}}$$

(12)

$$\nabla_{\mathbf{e^{r}}}\mathbf{L}=\frac{\partial\mathbf{L}}{\partial\mathbf{P}}\frac{\partial\mathbf{P}}{\partial\mathbf{H^{r}}}\frac{\partial\mathbf{H^{r}}}{\partial\mathbf{e}^{r}}$$

Here L is the training loss and P includes link prediction scores for vertices. The dimension of P and L is $1\times|V|$.

Although both $\mathsf{HDReason}$ and previous GNN model (schlichtkrull2018modeling, ; vashishth2019composition, ; shang2019end, ) use front-end encoder (graph model) and backend score function structure, the base hypervector matrix remains fixed in hyperdimensional computation. This leads to more efficient learning since HDC only updates the original embedding of vertex and relation.

Compared to conventional GNN models (schlichtkrull2018modeling, ; vashishth2019composition, ; shang2019end, ), HDC boasts enhanced memorization and information-retrieval proficiencies, allowing for the reconstruction of neighboring data for each vertex, as showcased in (poduval2022graphd, ). As a result, the memorization hypervector ($\textbf{M}^{\textbf{v}}$) of each vertex symbolically represents its neighboring vertex’s details. In contrast to the standard GNN model, $\mathsf{HDReason}$ offers superior interpretability and is more aligned with the principles of artificial general intelligence (AGI) as indicated in (latapie2021metamodel, ).

This section introduces the architecture of $\mathsf{HDReason}$, tailored for a CPU-FPGA heterogeneous computing platform, as illustrated in Figure 3. Our design is comprised of two parts: (1) the CPU component (host side) and (2) the FPGA component (kernel side), the latter one is where the original vertex and relation embeddings are managed and updated. The CPU-based scheduler is responsible for offloading vertex and relation data to the FPGA kernel. Within the FPGA kernel, the Encoder IP maps the vertex and relation embedding from the original space into the hyperspace and then performs the HDC memorization process. Every vertex’s resulting memorization hypervector is stored in the high bandwidth memory (HBM). The Score function IP then calculates and dispatches the reasoning outcome back to the host CPU. During the backpropagation phase, the bulk of gradient computations are offloaded to the Training IP.

The KG memorization includes two parts. The first part is an out-of-order (OoO) scheduler running on CPU trying to balance different vertices’ computations. The second part is a hardware acceleration IP conducting hypervector encoding and vertex neighbor aggregation.

After generating memory hypervectors ($\textbf{M}^{\textbf{v}}$), the next step is to perform reasoning tasks over specific triples. For each inference or training epoch, suppose the batch size is $|B|$. As is shown in Figure 6.(a), at the initial stages ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$a$}}\cr}}}$ and $\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$b$}}\cr}}}$), we load query vertex hypervector matrix ($\textbf{M}_{\textbf{q}}^{\textbf{v}}$) and query relation hypervector matrix ($\textbf{M}_{\textbf{q}}^{\textbf{r}}$) onto the on-chip buffer from FPGA storage (HBM). As our goal is to predict whether each triple exists or not, we pipeline the score function calculation process and only load one memory hypervector at a time ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$c$}}\cr}}}$), instead of loading the entire memory hypervector ($M^{v}$). In Figure 6, we assume that the Score function IP is evaluating vertex i over the query batch. Since the query batch contains $|B|$ subject and relation combinations, we replicate the loaded memory hypervector into $|B|$ on-chip buffers ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$d$}}\cr}}}$). Next, the query hypervectors ($M_{q}^{v}$ and $M_{q}^{r}$) and the memory hypervectors ($M^{v}_{i}$) are loaded into $|B|$ Score Engine units, with each unit calculating the corresponding prediction score for each query.

Figure 6(b) presents a single Score Engine IP’s microarchitecture design. Suppose the illustrated Score Engine unit in Figure 6(b) has an index of $j$, then it calculates $j^{th}$ vertex of the batch’s score over vertex $i$. The first step is to perform hypervector addition over subject HDV (query vertex HDV $M^{v}_{j}$) and relation HDV ($M^{r}_{j}$) to obtain the object HDV ($M^{o}_{j}$). The second step is to calculate the hypervector distance between $M^{o}_{j}$ and $M^{v}_{i}$. The intermediate prediction result($P^{I}_{j}$) is loaded into the L1 Norm IP for L1 normalization over the hyperspace, resulting in the normalized prediction ($N^{p}_{j}$) and its gradient with respect to the queried vertex and relation hypervector. Notice that these two gradients are actually equal in value since they are added together in the score function: $\frac{\partial N^{p}_{j}}{\partial M^{v}_{j}}=\frac{\partial N^{p}_{j}}{\partial M^{r}_{j}}$. As mentioned in section 4.2, to reduce the computation cost of back-propagation and avoid unnecessary data movement, inside the L1 Norm IP, we perform gradient computation of $\frac{\partial N^{p}_{j}}{\partial M^{v}_{j}}$ during the forward propagation, as shown in Figure 6(c) and Figure 6(d). Inside L1 Norm IP, there are $D$ Norm Units in total. Each of them extracts the absolute values and signs of elements of the input hypervector. The results are then fed into Tree Adder IP and a concatenate IP (Concat IP) to generate norm prediction $N^{p}_{j}$ and $\frac{\partial N^{p}_{j}}{\partial M^{v}_{j}}$ respectively.

Until now, after the execution of Score Engine IP, we have each batch member’s both forward path and backward path results. Next we pass the forward path result ($N^{p}_{j}$) and backward gradient result ($\frac{\partial N^{p}_{j}}{\partial M^{v}_{j}}$) to a Concat IP ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$e$}}\cr}}}$) and Tree Adder IP ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$f$}}\cr}}}$) respectively. To generate the final result, we need to concatenate each batch member’s L1 norm result ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$g$}}\cr}}}$) and accumulate all batch member’s gradient hypervectors ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$h$}}\cr}}}$). In the end, we will pass the forward path result ($N^{p}$) back to the CPU for further processing, such as the Sigmoid function ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$i$}}\cr}}}$). The gradient result ($\frac{\partial N^{p}_{j}}{\partial M^{v}_{j}}$) instead, will be stored inside the FPGA storage (such as HBM) for later training back-propagation usage ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$j$}}\cr}}}$).

Figure 7 presents the vertex embedding training computing flow. For the readability of the figure, we omit the detail of relation embedding training as it is similar to Figure 7. For illustration convenience, we rewrite the backpropagation equation:

(14)

$$\nabla_{\mathbf{e^{v}}}\mathbf{L}=\frac{\partial\mathbf{L}}{\partial\mathbf{N}}\frac{\partial\mathbf{N}}{\partial\mathbf{N^{p}}}\frac{\partial\mathbf{N^{p}}}{\partial\mathbf{M^{v}}}\frac{\partial\mathbf{M^{v}}}{\partial\mathbf{H^{v}}}\frac{\partial\mathbf{H^{v}}}{\partial\mathbf{e}^{v}}$$

We first rewrite the first two terms as:

(15)

$$\mathbf{\delta^{v}}=\frac{\partial\mathbf{L}}{\partial\mathbf{N}}\frac{\partial\mathbf{N}}{\partial\mathbf{N^{p}}}$$

We will calculate each batch member’s $\delta^{v}$ on CPU and concatenate them together. The generated new matrix is named $\mathbf{\delta}$ whose dimension is $|B|\times|V|$ where $|B|$ is the batch size and $|V|$ is the vertex size of the graph. Since graph size $|V|$ could be very large (over 10K), we cut $\delta$ into several chunks ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$A$}}\cr}}}$). The size of each chunk is $|B|\times T$. In Figure 7, we suppose the chunk index is i. So this chunk could be represented as $\delta[i*T:(i+1)*T-1]$. We also pipeline the computation between different chunks which means we will update T vertex embedding at each period.

Next we pass the $\delta[i*T:(i+1)*T-1]$ to the kernel FPGA ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$B$}}\cr}}}$). In section 4.2 and 4.3, we already compute the gradient result of $\frac{\partial N^{p}}{\partial M^{v}}$ and $\frac{\partial M^{v}}{\partial H^{v}}$ during the forward propagation and stored them inside the FPGA storage (such as HBM). Also the gradient of $\frac{\partial\mathbf{H^{v}}}{\partial\mathbf{e}^{v}}$ is equal to the transpose of base HDV ($\textbf{H}^{\textbf{B}}$). We only need to load precomputed matrices into the on-chip buffer ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$C$}}\cr}}}$, $\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$E$}}\cr}}}$, and $\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$G$}}\cr}}}$) and use two systolic arrays (SA) IP and one element-wise multiplication (MUL) IP to perform the actual computation ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$B$}}\cr}}}$, $\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$D$}}\cr}}}$, and $\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$F$}}\cr}}}$). After the computation of systolic array IP 2 (SA2), we will get T vertex’s gradient result $\nabla$L which will be passed back to the host CPU ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$H$}}\cr}}}$). The last step is updating the T vertex embedding model ($\mathbin{\vphantom{\circledast}\text{\ooalign{\smash{\raisebox{-2.58334pt}{\scalebox{2.3}{$\bullet$}}}\cr\smash{{\color[rgb]{1,1,1}\bf\footnotesize$I$}}\cr}}}$) which is also the last stage of the whole pipeline computation.

We benchmark $\mathsf{HDReason}$’s KGC accuracy over two standard datasets: FB15K-237 (toutanova2015observed, ) and WN18RR (dettmers2018convolutional, ) which are widely used for graph-related tasks. For $\mathsf{HDReason}$’s hardware acceleration performance, we also added two large-scale datasets: WN18 (miller1995wordnet, ) and YAGO3-10 (mahdisoltani2014yago3, ). Table 3 lists the properties of these four KG datasets. Our platform comprises an Intel i9-12900KF CPU as the host, along with a Xilinx Alveo U280 (or U50 depending on the configuration) card to implement the FPGA kernels. We implement the kernel accelerator using SystemVerilog and synthesize it using Xilinx Vivado. The communication between the host CPU and kernel FPGA is supported by the Xilinx Vitis platform via PCIe direct memory access (PCIe DMA) (kathail2020xilinx, ).

Figure 8.(a) and Figure 8.(b) present the $\mathsf{HDReason}$’s reasoning accuracy over FB15K-237 and WN18RR. For double direction reasoning, we compare our HDC-based model with other state-of-the-art graph models including TransE (bordes2013translating, ), R-GCN (schlichtkrull2018modeling, ), SACN (shang2019end, ), and Comp-GCN (vashishth2019composition, ). For single direction reasoning, we compare our HDC-based model with other state-of-the-art reinforcement learning (RL) models including MINERVA (das2017go, ), C-MINERVA (stoica2020contextual, ), R2D2 (hildebrandt2020reasoning, ), RARL (wu2021findings, ), and ADRL (wang2020adrl, ). Table 4 summarizes critical model parameters including original vertex embedding dimension d and final embedding dimension D (after encoding and aggregation). Figure 8.(a) shows that $\mathsf{HDReason}$ achieves relatively high reasoning accuracy compared to state-of-the-art GCNs such as CompGCN and SACN. Besides, $\mathsf{HDReason}$ only needs to train vertex and relation embedding, which simplifies the model training computation complexity.

Figure 9.(a) illustrates the HDR model’s high robustness during KGC tasks. After encoding and memorization, we attempt to drop some hypervector dimensions before entering the score function. HDC model achieves significant robustness when dropping dimensions with low entropy. This attribute is also reported by previous HDC work (zou2021scalable, ). This characteristic allows the HDC model to be deployed on resource-constrained platforms by maintaining only a portion of the model’s hypervector weight. On the other hand, Figure 9.(b) demonstrates that, when compared to a traditional GNN-based model (specifically SACN), HDR supports aggressive low-bit quantization. Here we try to quantize HDR and SACN into different fixed-point numbers using QPyTorch (zhang2019qpytorch, ). When quantizing the HDR model from floating-point numbers to 4-bit fixed-point numbers, it experiences only a 5% drop in accuracy, whereas the GNN-based model experiences a 45% accuracy drop. This indicates that HDR is well-suited for low-bit precision computation platforms, such as FPGA.

Table 5 presents the FPGA resource utilization of $\mathsf{HDReason}$. The vertex initial embedding dimension (d) is 96 and the encoded hypervector dimension (D) is 256. The training batch size is 128 and pipeline chunk size T is 32. We use 8 HBM pseudo channels (PCs) and AXI with 256-bit data width: 4 PCs are used to store vertex ($\textbf{H}^{\textbf{v}}$) and memorization hypervectors ($\textbf{M}^{\textbf{v}}$); while another 4 PCs hold gradients including $\frac{\partial N^{p}}{\partial M^{v}}$ and $\frac{\partial M^{v}}{\partial H^{v}}$. We save relation hypervectors ($\textbf{H}^{\textbf{r}}$) in the on-chip UltraRAM. The ”Other” section in Table 5 includes all the Xilinx IP resource usage such as AXI Interconnect IP and PCIe DMA (feist2012vivado, ).

Table 6 tabulates the execution time, energy, and memory usage of $\mathsf{HDReason}$ on FPGA and GPU. Our GPU implementation uses PyTorch (paszke2019pytorch, ) and PyG (fey2019fast, ). For a fair comparison, we specify both FPGA and GPU’s training batch size as 128 except YAGO3-10 on GPU. We collect the FPGA power using Xilinx Power Estimator (XPE) and GPU power using NVIDIA Management Library (NVML). When testing $\mathsf{HDReason}$ on RTX 3090 using the YAGO3-10 dataset, due to limited GPU memory (24 GB) and large-scale graph vertex size (over 120K vertices), we decrease the training batch size from 128 to 32. As for speedup, our accelerator shows, on average, over 9$\times$ improvement compared to GPU even Alveo U50 adopts less advanced fabrication technology compared to RTX 3090 (14 nm FinFET on FPGA vs 8 nm TSMC on GPU). Our FPGA-based accelerator shows notable advantages over GPU in terms of speed, energy, and memory efficiency. This large advantage comes from our hardware optimization shown in Figure 8.(c), including reusing encoded hypervectors, a density-aware scheduler to balance the computation, and computing backward gradients in the forward path. For energy efficiency, our design shows over 95$\times$ improvement compared to RTX 3090. For memory efficiency, as we cut the large vertex embedding gradient computation into small chunks and pipeline computation between different chunks, we avoid using unnecessary intermediate memory on the FPGA. Therefore, our FPGA accelerator supports a high batch size on large-scale graph datasets (like YAGO3-10).

Figure 8.(d) shows the execution time breakdown during each single batch training. Here, CPU time includes data communication time between host CPU and kernel FPGA, CPU gradient computation (equation 15), and vertex and relation embedding update. Since we compute partial backward gradient in the forward path, in Figure 8.(d), actual training computation (Training) only takes a very small part of the overall training process. The memorization (Mem) execution time instead takes over 50% total execution time. In section 5.5, we will discuss how to decrease the memorization execution time by tuning the accelerator’s parameters.

In this section, we discuss the effects of vertex replacement policy and on-chip storage usage over memorization speed. Compared to previous GCN acceleration works, our accelerator reuses encoded hypervectors such that the main computation overhead switches from matrix multiplication to the data transfer between FPGA and HBM. In Figure 10, we show the changes in memorization speed and FPGA HBM data communication with different policies and memory usage. In Figure 10, those UltraRAMs are only used to store vertex hypervectors ($\textbf{H}^{\textbf{v}}$). Due to the limited on-chip storage, the FPGA will need to replace hypervectors currently stored on-chip with what the incoming ones needed for performing aggregation operations. This process is like a cache miss in the CPU memory hierarchy. Those access misses will bring unnecessary data communication between FPGA and HBM. As is shown in Figure 10, when on-chip UltraRAM usage is small (such as 64 UltraRAM), then data transfer between FPGA and HBM is large. As expected, with the increased usage of on-chip UltraRAM, both the memorization time and FPGA-HBM communication time decrease. This trend is especially noticeable for large-scale datasets like YAGO3-10 (Figure 10.(d)).

In Figure 10, we present the influence of different replacement policies. We tried three different replacement policies including random replacement (Random), last frequently use (LFU), and last recent use (LRU). On average, LFU achieves the best performance and 8% better than Random. To further extend our findings, future work could involve designing specific replacement policies targeting KG datasets.

Figure 11 presents the speedup and energy efficiency comparison between different models on different hardware platforms. Specifically, we show 4 different graph models on 8 different hardware platforms. For 8 hardware platforms, we try 2 different CPUs (Intel i9-12900 KF and AMD Ryzen 5955 WX), 3 different GPUs (NVIDIA RTX 3090, RTX 4090, and A100), and 5 different FPGAs (Xilinx Alveo U50, U200, U250, U280, and Kintex7 KC705). On the GPU platforms, we choose to use PyG to optimize models’ execution speed and GPU memory usage.

For FPGA platforms, the configuration of our design on Alveo U50 is the same as Table 5. For Alveo U280, with more resources, we increase the HBM channels from 8 PCs to 16 PCs, set the AXI data width as 512, increase the number of Memorization Computing IP ($\textbf{N}_{\textbf{c}}$) from 16 to 32, training chunk size (T) from 32 to 64, and use 256 UltraRAMs to store vertex hypervectors on-chip. We approximate the HDR’s performance on existing state-of-the-art HDC FPGA accelerator (LookHD) (imani2021revisiting, ). We also approximate the performance of 4 different models on Alveo U200 and U250 based on state-of-the-art GNN training FPGA acceleration works (GraphACT (zeng2020graphact, ) and HP-GNN (lin2022hp, )).

To ensure a fair comparison, we compare the small FPGA design (Alveo U50) with GraphACT which is based on Alveo U200, and compare the large FPGA design (Alveo U280) with HP-GNN which is based on Alveo U250. Compared to GraphACT (Alveo U200), $\mathsf{HDReason}$ on Alveo U50 shows on average 9$\times$ speedup and 10$\times$ energy efficiency improvement. When comparing with HP-GNN (Alveo U250), $\mathsf{HDReason}$ on Alveo U280 shows on average 3.5$\times$ speedup and 4.6$\times$ energy efficiency improvement. When comparing with RTX 4090, our Alveo U280 accelerator achieves on average 10.6$\times$ speedup and 65$\times$ more energy efficiency. When conducting cross models and platforms comparison, $\mathsf{HDReason}$ on Xilinx Alveo U280 (and Xilinx Alveo U50) provides on average 4.2$\times$ (and 8.3$\times$) speedup and 3.4$\times$ (and 9.1$\times$) energy efficiency improvement with similar accuracy as compared to the GCN training platform HP-GNN (and GraphACT).