Across Time and Space: Senju’s Approach for Scaling Iterative Stencil Loop Accelerators on Single and Multiple FPGAs

Because of their iterative nature, ISLs usually represent the bottleneck of a given application, making them good candidates for hardware acceleration. In addition, the memory transfers frequently bound the performance of these algorithms, as they have to access multiple data every time to perform the stencil calculation. Consequently, a critical part of ISL acceleration is the design of an efficient memory sub-system that guarantees continuous data processing. In this scenario, FPGAs have an advantage over other architectures [2, 9, 21, 36, 52] as they allow building a custom sub-system for the target computation. For instance, Richter et al. [43] proposed a flexible RTL template for 3D stencils employing local yet large on-chip buffers to store input data and minimize off-chip transfers. Conversely, the ISL automatic flow by Nacci et al. [38] addresses memory boundness by creating a cone-shaped structure where subsequent tiles process multiple merged iterations concurrently. However, such an approach requires redundant computations among neighboring tiles and causes on-chip memory contention, negatively affecting performance. The OpenCL-based framework by Wang et al. [58] employs a similar approach and avoids redundant computations through pipes between tiles. Nonetheless, this solution does not solve the on-chip memory contention, reducing the achievable performance as the number of merged iterations grows.

Literature FPGA-based architectures mainly rely on spatial and temporal parallelism to boost ISL performance. The former refers to the amount of parallel input data processed per stencil module, while the latter orthogonally refers to the deepness of the stencil module pipeline. For instance, Figure 2(a) illustrates a stencil module with a vectorized packet of s words, while Figure 2(b) represents a queue of t stencil module replicas performing t timesteps with a single input data word. For instance, Kobori and Maruyama leveraged temporal parallelism to accelerate ISLs for cellular automata [34]. On the other hand, the SODA framework [8] builds a dataflow microarchitecture that exploits both parallelism types and minimizes external memory transfers and on-chip memory usage. Furthermore, various studies exploit multi-FPGA approaches to increase performance through temporal parallelism. For example, Sano et al. [45] described a flexible, runtime-configurable multi-FPGA architecture for ISLs, showcasing linear performance scalability on a cluster of nine FPGAs. Likewise, Waidyasooriya et al. [56] proposed a multi-FPGA architecture employing both temporal and spatial parallelism. They evaluated the performance of their two-FPGA system, which builds upon a previous one by the same authors [57], on 2D and 3D benchmarks using skewed and non-skewed grids. Finally, SASA [54] is a framework that exploits HBM banks available on recent FPGAs to implement hybrid spatial and temporal parallelism. To this end, the authors designed two ways to reduce the redundant computations derived from their hybrid approach and evaluated different parallelism combinations.

SST Microarchitecture. The SST microarchitecture [6] defines an exciting approach for ISL acceleration on FPGA. Its memory sub-system employs a non-uniform memory partitioning [11] to achieve on-chip buffering and supply parallel access to stencil points at each clock cycle. Specifically, the sub-system relies on a series of hardware modules, one per stencil point, that filter the incoming data according to internal counters and dispatch them to the following ones through on-chip FIFO buffers. Eventually, the data reach the module that implements the stencil computation. Finally, the last module reconstructs the spatial coherency of the input, enabling SST queuing [6]. SSTs went through various incarnations [6, 40, 42], and each research work added features to improve the performance of both single- and multi-FPGA implementations. For instance, the latest version of SST microarchitecture described by Reggiani et al. [42] introduced spatial parallelism, even though they limited it to four parallel stencil computations. Nonetheless, we believe there is still room for significant improvements at different levels (e.g., parallelism, microarchitecture, and design automation) for such a relevant methodology for ISL accelerators.

Our Proposal. Given this context, Senju represents the latest SST incarnation, presenting new features and enhancements. At the microarchitecture level (Section 3), we introduce various improvements to reduce resource usage and increase performance. In particular, we implement a hybrid FIFO approach combining single- and packed-data FIFOs (Section 3.2), which differs from literature solutions employing single-data (e.g., SODA [8] and Reggiani et al. [42]) or packed-data FIFOs only (e.g., SASA [54]). Then, compared to previous SST solutions [6, 40, 42], we optimize and reduce the internal logic usage and increase the supported spatial parallelism. Moreover, unlike most literature studies focusing on frameworks for the single-FPGA scenario [8, 54] or manually implementing multi-FPGA designs [42, 45], we offer a comprehensive design automation framework for single-/multi-FPGA systems featuring a performance model and multiple configuration parameters (Section 4). Finally, we expand the range of experimental evaluations and literature comparisons to provide a deeper analysis of our results through unexplored metrics and methodologies (Section 5).

An SST is a hardware module that performs a timestep of an ISL operating on an input grid of points. Specifically, this approach relies on a dataflow computational model that decouples the management of each stencil point, enabling parallel access per clock cycle. To this end, an SST comprises four fundamental fully-pipelined components: Filters, FIFO buffers, Processing Element (PE), and Multiplexer (MUX). Figure 3 shows a comprehensive example of an SST module for a 3-point Jacobi 1D stencil, which includes the computation pseudocode (Listing 1), the content of input and output arrays (Figure 3(a)), and the SST microarchitecture (Figure 3(b)). We will use Figure 3 as a reference throughout this section to ease the comprehension of the concept behind the SST design.

SST Filter. Each SST contains as many Filters as the number of points the target stencil processes. Given a mapping between a specific stencil point and a Filter, the role of this component is to dispatch input data (coming from another hardware module or a previous Filter) to the following components. Specifically, a Filter communicates with the next one (except for the last Filter) and the PE. In addition, a specific Filter, usually the central one, also connects to the MUX and transmits the input boundary cells. While the communication between Filters always happens, the one with other components relies on filtering conditions ruled by internal input counters, whose number equals the input grid dimensionality. Thus, as soon as enough data have passed through the Filter, the counters activate the conditions and enable data dispatching to the PE (or MUX). Please note that these conditions depend on the size of each input dimension and the stencil point the Filter maps onto. In Figure 3(b), the 3-point Jacobi 1D requires three Filters, each one assigned to a given point (e.g., F0 to in[i+1]) and enforcing different conditions (reported on the connections between Filters and PE/MUX) based on a single counter per Filter (e.g., c0_F0 for F0) since this stencil is 1D.

SST FIFO. A peculiar feature of the SST microarchitecture is the memory sub-system, which implements a non-uniform partitioning scheme, where on-chip FIFOs placed between internal SST components buffer uneven input portions and enable data reuse. This scheme derives from internal synchronization mechanisms that some components enforce, preventing a seamless data flow. Specifically, PE and MUX require concurrent and conditioned access to their inputs, respectively, making data accumulate inside such FIFOs in the meantime. Conversely, the smooth data transfer among Filters demands relatively small buffers. Hence, each FIFO depth has to be carefully specified to avoid stalls/deadlocks or resource waste. In Figure 3(b), we indicate the content of FIFOs to ease the understanding of the data flow. Due to the 1D nature of this stencil and the limited number of points, small FIFOs are enough to avoid stalls and guarantee a dataflow execution.

SST PE. The PE is the component in charge of performing the stencil computation. Specifically, the PE features as many inputs as the number of stencil points, which come from the Filters and wait inside the FIFOs. Once all the points are available, the PE concurrently accesses/reads them from the FIFOs, executes the stencil calculation, and forwards the result to the MUX. In Figure 3(b), we depict how the PE receives the input stencil points from the three Filters, performs the operation described in lines 3–4 of Listing 1, and finally transmits the output to the following component.

SST MUX. The last component of an SST is the MUX, which receives inputs from both the central Filter and the PE and outputs spatially coherent data; in other words, the data position in the output array is the same as the input one. The MUX obtains this result through selective read access to the input FIFOs based on internal conditions and input counters (one per input dimensionality, resembling the Filter’s mechanism. Consequently, when the MUX reads data from one component (e.g., the central Filter) through the FIFO, further data from the other (e.g., the PE) accumulate in the corresponding interposed FIFO. In Figure 3(b), we show that the MUX produces spatially coherent results (i.e., in[0], out[1], out[2], and so on) out of the inputs from PE and F1 and based on the internal conditions and counter c0_MUX.

Scaling Input Grid Dimensionality. Although we exploited a 1D stencil as a reference (Figure 3), the SST microarchitecture can be easily generalized to a higher input dimensionality. Indeed, the dimensionality directly affects the internal conditions and counters of Filters and MUX; conversely, other SST features (e.g., the number of Filters and FIFO sizes) still depend on the stencil shape. For instance, a 5-point cross-shaped 2D stencil, such as the one in Figure 1, would require five Filters (one per point), a five-input PE, and FIFO buffers large enough to potentially accommodate rows of the input grid since the stencil spans three rows. Moreover, due to the bi-dimensional input, each Filter and the MUX implement more complex conditions based on two internal counters.

Exploiting Temporal Parallelism. The role of the MUX is paramount in the big picture of ISL acceleration through SSTs as it enables exploiting temporal parallelism, that is, serially connecting various SST-based modules and building a deep pipeline that accelerates multiple timesteps in a single sweep. To this end, each SST module implements a streaming interface (Intel’s Avalon in our case [25]) featuring valid, ready, and data signals. This way, the overall pipeline dynamically manages the data flow through a backpressure mechanism. In addition, to provide data to this pipeline, we need two additional interface modules: the input interface module reads data from a source and sends them to the first SST, and the output one forwards the results of the last SST to a given destination. In our context, the source/destination may be the off-chip memory or another FPGA. In the former case, the communication happens via a memory-mapped interface with the off-chip memory and a streaming one with the stencil pipeline [25]. Conversely, the latter case employs streaming interfaces for input and output, requiring two additional signals, startofpacket and endofpacket, when dealing with other FPGAs.

The microarchitecture we defined so far executes a single stencil computation per clock cycle and supports temporal parallelism by queuing multiple SST modules. On the other hand, spatial parallelism represents another orthogonal way to improve performance by simultaneously computing the same stencil on multiple data. To this end, spatial parallelism requires feeding the SST-based module with more than one single input word per clock cycle. Thus, the interface modules must read/write vectorized/packed data, that is, packets containing multiple words (for instance, a 512-bit packed input includes 16 32-bit floating-point words), from/to sources/destinations, increasing the overall bandwidth utilization until saturation. Similarly, each SST component must correctly handle, dispatch, process or unpack such packed data through a revised internal structure (e.g., conditions and counters). Overall, spatial parallelism seriously impacts each SST’s resource usage.

Given these motivations, we now describe how we adapted and optimized the basic SST module for more spatial parallelism and the various enhancements we implemented to reduce resource footprint. To this end, unlike Section 3.1, we take as an example a 5-point Jacobi 2D stencil with spatial parallelism = 16 because its greater complexity than the 1D case enables us to highlight our optimizations’ benefits better; nonetheless, our considerations also hold for other stencils and input dimensionalities. We reimplemented the SST architecture by Reggiani et al. [42] through the Intel HLS Compiler and employed it as a baseline to assess how our incremental optimizations quantitatively affect resource usage, frequency, and the theoretical maximum temporal parallelism, as reported in Table 1. Specifically, even though the Adaptive Logic Module (ALM) amount already includes Memory Logic Array Blocks (MLABs) (usually a quarter of ALMs [32]), we include both resources due to the latter’s importance in implementing FIFOs. Indeed, when a FIFO depth is relatively small, MLABs store its content. Conversely, large FIFOs employ RAMs, allowing MLABs to serve as ALMs. Resource and frequency values come from hardware syntheses of the different HLS designs performed by Intel’s Quartus [29] targeting an Intel Stratix 10 GX FPGA [30]. Finally, temporal parallelism depends on the critical resource of each design (in red).

Packing Input Filter FIFOs. The first significant design choice involves the management of the packed input data when arriving at the SST-based module. Indeed, even though the Filters generally behave similarly to their counterparts we previously described, their data dispatching mechanism differs. For instance, Reggiani et al.’s SST architecture [42] always works with unpacked data and employs separate single-data FIFOs (i.e., a FIFO containing unpacked data) to transmit single words among internal components. In other terms, if the input packet contains four words, each Filter uses four FIFOs to forward such values to the next Filter. However, there is no need to unpack and separately transfer the data from one Filter to another, as no filtering operation occurs on the single words for this communication. Thus, our architecture features a single packed-data FIFO (i.e., a FIFO containing packed data) to distribute packed data among Filters. In this way, we decrease the usage of ALMs and Registers thanks to the reduction/removal of redundant control within Filters. Finally, RAM usage slightly reduces, but MLAB demand grows because of the internal policies the Intel HLS Compiler enforces for FIFO mapping, even if it remains below 1.5% (Table 1, version 2).

Single PE and MUX. Moving to other SST components, as the spatial parallelism increases, one way to deal with it consists of instantiating N PEs/MUXs to execute N operations concurrently [42]. However, this design choice implies replicating common logic among the components. Thus, our solution implements a single PE and MUX. The former component receives data from all the Filters, computes the same stencil N times on different inputs in parallel, and sends the output to the MUX. The latter reads the input words coming from the PE and the central Filter FIFOs, combines them into a packet, and then transmits it to the following hardware module. This way, the HLS compiler can optimize and reuse common logic and save resources (Table 1, version 3).

Hybrid FIFO Approach. Moving to the communication between Filters and PE/MUX, a Filter cannot dispatch a packed input data to these components despite its revised design. Indeed, although the spatial parallelism growth enables relaxing the filtering conditions, there are still some problematic/corner cases (e.g., close to the borders of the input grid) that the components have to handle separately to avoid computation errors. For instance, as explained in Section 3.1, each Filter dispatches a specific subset of stencil points to the PE and MUX. However, when we move from a scalar case to a vectorized one, packed data may comprise points outside this subset. Besides, once all the inputs are available, the PE may immediately process only a portion of their points and have to buffer and reorganize the remaining ones while waiting for the availability of the following inputs, complicating the internal logic. In the SST scenario, Reggiani et al. [42] addressed this issue through separate single-data FIFOs, one per word inside the packed data; conversely, literature studies based on non-SST designs use single- or packed-data FIFOs [8, 54]. Instead, we propose a novel hybrid FIFO approach to manage the data transfer between Filters and PE, central Filter and MUX, and PE and MUX. Specifically, a hybrid FIFO is a logical combination of single- and packed-data FIFOs; given packed input data, a Filter uses single-data FIFOs to dispatch aforementioned corner-case words and a packed-data FIFO to transfer the remaining ones, whose number increases with the spatial parallelism. This way, the next components (i.e., PE and MUX) can decouple the input management, independently process the data from both types of FIFO, and, in the case of the PE, dispatch the results to the MUX. Of course, each component enforces new internal conditions to manage such cases properly but reutilizes the previous counters. Besides, even the PE now employs counters (unnecessary in the scalar case). This hybrid approach notably reduces resource usage, especially ALMs, making RAMs the critical one (Table 1, version 4).

Stall-enabled Cluster. Finally, the last optimization directly derives from how the Intel HLS Compiler works. Indeed, this compiler generally packs related operations into clusters, mainly stall-free [24]. The advantage of this cluster type is that the target operations execute without stalls. Besides, the cluster includes a FIFO at its end to hold the results if it is stalled. This way, the resulting design may achieve a higher running frequency and throughput at the cost of area and latency. However, if we intend to reduce resource usage further and can tolerate a frequency lowering, we can employ the hls_use_stall_enable attribute to disable this feature and force the compiler to generate stall-enabled clusters that do not contain FIFOs. Hence, we obtain a relevant reduction for most resources but MLABs, which become the critical one (Table 1, version 5).

Summary. We enhanced SST design through microarchitectural optimizations to better support spatial parallelism and reduce resource usage, paving the way to more extensive temporal parallelism. Indeed, in our example, we increased the theoretical temporal parallelism from 15 to 55; please note that, due to the relevant impact of the last proposed optimization, Senju does not disable stall-free clusters by default, offering the user the option to do that. Figure 4 depicts the optimized SST microarchitecture for the 3-point Jacobi 1D example from Section 3.1, showcasing the internal changes to support the vectorized scenario (i.e., spatial parallelism > 1).

We now focus on how a stencil accelerator, which features the SST pipeline and interface modules, integrates within an overall system. Currently, Senju targets the PCIe-based Programmable Acceleration Cards (PACs) [27] by Intel. The FPGA on such boards is usually organized in two regions: a static one and a partially reconfigurable one. The former contains the FPGA Interface Manager (FIM) [26], which manages the interaction and communication with the various off-chip components available on the board (e.g., PCIe, DDR banks, and QSFP28). The latter contains the hardware design of the application to accelerate, usually called Accelerator Functional Unit (AFU). As the FIM does not implement a direct connection between the PCIe and other components on the board, for example, the DDR banks, we rely on a shell within the AFU region for such a purpose. Such a shell [46, 55] features multiple components to interface an accelerator with the off-chip memory or other FPGAs. For instance, it includes a DMA controller to write data coming from the host machine through PCIe to the off-chip memory and vice versa. Similarly, the shell leverages Flow Controllers (FCs) [37] and serial transceiver modules, namely Serial Lite III (SL3) IPs, to implement point-to-point inter-FPGA communication through the PAC network ports, which we will exploit to build multi-FPGA systems based on ring or chain topologies (see Section 4.4). In particular, the FC implements backpressure functionality to prevent buffer overflows, which, on the other hand, causes a slight throughput degradation [37]. Finally, the shell contains a \(4\times 4\) memory-mapped crossbar to route the internal data traffic based on a runtime configuration (for instance, from the off-chip memory to the network, passing through the DMA). Given this infrastructure, the stencil accelerator can communicate with the off-chip memory or another FPGA according to the topology Senju implements, deriving from the user’s specification (see Section 4.4). Figure 5 shows the overall system, including the shell and the stencil accelerator in the AFU.

Senju accepts as input a JSON file that describes the stencil computation to implement and includes additional information to guide the design process. Specifically, the list of expected/optional fields in the JSON file is the following:

After parsing such a file and checking its correctness, Senju invokes its internal modules and builds the design according to the user’s specification. To exemplify, Listing 1 shows a possible JSON file structure. In this case, we require Senju to implement a 5-point Jacobi 2D stencil processing \(1024\times 1024\) floating-point data. Besides, the input description file indicates the spatial parallelism but not the optional temporal parallelism; hence, Senju will automatically calculate it based on the resource usage of the requested stencil. Then, the JSON file defines the target board, skips the design synthesis, and requires performance estimation through Senju’s model. Finally, the input file specifies the remaining fields, namely frequency, the multi-FPGA setup (i.e., number of FPGAs and topology), and the usage of cluster optimization.

Once Senju has received and parsed the input JSON description, the first module it invokes is the Stencil Generator ①. Such a module is in charge of generating the hardware architecture for the target stencil according to the user’s specifications. Specifically, it requires the following JSON fields: stencil, input size, spatial parallelism, topology, and cluster optimization. On the other hand, the outcome of this module is set of C++ codes suitable for the Intel HLS Compiler.

If the user does not specify the temporal parallelism they want to implement for the target stencil computation, Senju invokes the Stencil Counter ② module. Such a module takes the code generated by the previous one and employs the HLS compiler to produce the Register Transfer Level (RTL) of the stencil component. Then, it utilizes the synthesizer to compile RTL and extracts the resource usage of a whole SST module, which includes the internal components and FIFOs. At that point, the Stencil Counter module can calculate the maximum temporal parallelism. To this end, the module first computes how many usable resources the FPGA provides to instantiate the stencil pipeline:where R is the set of FPGA resources (ALMs, Registers, RAMs, DSPs, and MLABs), \(FPGA_{r}^{u}\) and \(FPGA_{r}\) indicate the usable and available amount of resource r on the target FPGA, respectively, and \(system_{r}\) is the quantity of resource r required by other components on the FPGA (e.g., FIM, shell, and input/output interface modules). Then, the module identifies the critical resource as follows:where cr is the critical resource, \(stencil_{r}\) is the demand of resource r by the stencil component, and budget is the resource budget the module employs to balance the temporal parallelism and the design routability. Consequently, temporal parallelism is:where the module removes one more stencil to relax design constraints further.

Finally, if the temporal parallelism is equal to or greater than a certain threshold (100 in our case), we increase the impact of the critical resource for routability reasons:

The FPGA Design Generator ③ module builds the overall hardware design to synthesize and automates all the steps toward the bitstream generation. Specifically, this module invokes the HLS compiler to produce the RTL of the stencil component (if the Stencil Counter module has not already done it) and the required interface modules. Then, the module builds one or more systems according to the number of FPGAs and topology that the user requested. In particular, we have two possible scenarios: single- or multi-FPGA designs. In the former case, the module instantiates the shell, the interface components to interact with the off-chip memory, and multiple stencils, whose number depends on temporal parallelism. Next, the module properly connects such components to build a deep pipeline of stencils that reads/writes data from/to the off-chip memory. Besides, when more than one DDR bank is available, the module uses separate banks to avoid congestion. In the latter case, the module has to potentially implement three designs (first, internal, and last) based on the chosen topology and the number of FPGAs. In the first FPGA of the chain topology, the module connects the input of the stencil accelerator (i.e., the stencil pipeline and interface modules) to the off-chip memory and the output to the 4 \(\times\) 4 crossbar to exploit the link with the SL3 IPs. Then, if the number of requested FPGAs > 2, the module implements a design for the internal FPGAs, which connects the input and output of the stencil accelerator to the 4 \(\times\) 4 crossbar to process the data coming from the previous FPGA and send the outcome to the next one. Finally, in the last FPGA design, the module connects the input of the accelerator to the crossbar and the output to the off-chip memory to store the results of the overall computation. Conversely, in the ring topology, the first and last FPGAs coincide. Thus, there is always at least one internal FPGA and the first FPGA also contains the logic to close the ring, that is, interface modules to read the data from the SL3 IP and store them in the DDR. Figure 8 shows an example of the two topologies, highlighting components involved in data transfer and processing.

Once the design process is over, the module starts the synthesis (using the default target frequency or the one suggested by the user) and waits for the bitstream generation. In particular, in the case of multi-FPGA systems, the module produces up to three bitstreams, that is, for the first, internal (if any), and last (if it does not coincide with the first) FPGAs.

The user can optionally ask Senju to estimate the Giga Operations per Second (GOPS) of the resulting single- or multi-FPGA design through the Performance Model ④ module. Such a simple yer effective module takes the input size and byte width, the untouched border of each input dimension (i.e., the boundary cells that the stencil does not modify), the target stencil operations, the target frequency, spatial and temporal parallelism, the off-chip memory and network bandwidths. Some of these values come from the user’s input description (e.g., input size and spatial parallelism), while others derive from the previous modules (e.g., boundary cells and stencil operations).

This module first calculates how many Giga Operations (GOPs) the stencil pipeline performs:where tp is the temporal parallelism, D is the input dimensionality, and \(input\;size_d\) and \(border_d\) are the size and untouched border of the input on dimension d, respectively, which define the number of cells to process. For instance, the 5-point Jacobi 2D stencil in Listing 1 performs five operations (the Intel HLS Compiler automatically evaluates the result of 1/5), and the \(border_d\) variable is 2 for both dimensions since the stencil does not affect the first and last rows and columns of the input grid, respectively (see Figures 1 and 7).

At this point, the module estimates the hardware design execution time. Specifically, it first calculates the total transferred data in gigabytes, which is twice the overall input size (including the border) because the amount of read/written data from/to the off-chip memory is the same:where \(data\;type\;bytes\) is the number of bytes of the input data type (e.g., 4 bytes for 32-bit floating-point data). Next, the module computes the execution time required to transfer such data from one DDR4 bank to another (or the same) through the FPGA, considering the bank bandwidth when the target frequency is f and spatial parallelism is sp (i.e., the number of words per packet): \(DDR\;bandwidth^{sp}_{f}\) depends on how many DDR4 banks we employ; as mentioned in Section 4.4, Senju leverages two banks, if available, to avoid reading/writing from/to the same bank, causing congestion and reducing performance. Thus, we extracted bandwidth values using one or two DDR4 banks by running multiple benchmarks with an increasing input size and spatial parallelism at different frequencies. We built a lookup table from these values, which the module consults according to the user’s input description. If the input size misses, the module interpolates the bandwidth, and, finally, scales it based on the target frequency the user specified, as long as this frequency does not oversaturate the bandwidth. In the case of a multi-FPGA design, the module employs the minimum between the network and off-chip memory bandwidths ( \(SL3\;bandwidth^{sp}_{f}\) and \(DDR\;bandwidth^{sp}_{f}\) , respectively), as it acts as the bottleneck:Here, to be compliant with \(transferred\;GB\) (which is twice \(input\;GB\) , but the amount of data passing through a network link is \(input\;GB\) ), we double \(SL3\;bandwidth^{sp}\) . Moreover, likewise the off-chip memory, we benchmarked the SL3 connection in various configurations (e.g., spatial parallelism and used DDR4 banks) and built another lookup table from the collected results.

In the last step, the model adjusts the transfer time by adding the delay introduced by each stencil in the accelerator, which depends on the temporal parallelism tp and estimates the GOPS:Since data from the central Filter and PE converge to the MUX, this component introduces an unavoidable delay in the SST module before producing results at a steady state due to the spatial coherency it must guarantee. Indeed, after outputting the initial border from the central Filter, the MUX waits for multiple clock cycles before receiving data from the PE. Specifically, as we empirically observed, the \(MUX\;delay^{sp}_{f}\) variable depends on such cycles, the MUX latency, the spatial parallelism sp, and the frequency f.

The first part of our evaluation focuses on analyzing the scalability of the designs generated by Senju as we increase both spatial and temporal parallelism. In particular, we are interested in measuring each design’s performance, energy efficiency, and resource usage scaling. Thus, on the one hand, we augment the temporal parallelism by ten stencils at a time until we reach the limit suggested by Senju; on the other, for each of the previous designs, we manually select five configurations of spatial parallelism (1, 2, 4, 8, and 16) to exploit the off-chip memory and network incrementally. Furthermore, for space reasons and to guarantee readability, we only report the Jacobi 1D, 2D, and 3D results scaling spatial parallelism values and employing up to 50 stencils for single-FPGA designs. On the other hand, for the same reasons, we focus on Jacobi 1D experiments with up to 100 stencils for multi-FPGA designs. Nonetheless, the same considerations also hold for greater temporal parallelism and other benchmarks (e.g., Heat 1D, 2D, and 3D). Finally, the input size is the same for these three benchmarks, namely \(2^{28}\) 32-bit floating-point data, organized as \(32768 \times 8192\) and \(32768 \times 128 \times 64\) grids for the 2D and 3D cases, respectively. This value guarantees a proper trade-off between off-chip memory/network bandwidth saturation and temporal parallelism.

This section first compares our designs against the ones available in the literature. Then, it assesses the quality of our results in opposition to another work based on the SST microarchitecture.

So far, we have discussed and analyzed the performance reached by Senju and other state-of-the-art studies in GFLOPS, GFLOPS/W, and bandwidth-normalized results. Indeed, the literature about ISL accelerators mainly concentrates on GFLOPS and GFLOPS/W metrics to assess the goodness of a given implementation, and, for sure, this approach provides valuable insights. Nonetheless, other methodologies are viable, as we proved when considering off-chip memory/network bandwidth. Similarly, we believe that obtaining the highest performance does not always imply the best results apriori, mainly when we target real applications that could use this kind of acceleration. In this scenario, the overall execution time depends not only on the accelerator performance but also on the number of times we invoke it. In particular, such a number hinges on how many iterations the application needs to converge. Thus, a proper balance between spatial and temporal parallelism is critical to reducing the overall application execution time. Indeed, high spatial parallelism lowers the accelerator execution time but demands more FPGA resources, diminishing temporal parallelism.

Let us consider an applicative scenario in line with the setup, we used for the previous experiments. In such a scenario, a real application could work like this in first approximation: (1) get the input data to process; (2) send them to the board via PCIe; (3) invoke the stencil accelerator; (4) read the data back from the board; (5) check the convergence; (6) if the check is successful, end the computation; otherwise, go back to step (2). In other words, if we exclude the first step, we can model the application execution time with the following formula:where \(tp_i\) and \(sp_j\) are temporal and spatial parallelism values, respectively, K is the transfer time between the host and the board, \(T_{tp_i, sp_j}\) is the stencil accelerator execution time when using a specific combination of \(tp_i\) and \(sp_j\) , C is the convergence check time, and It is the maximum number of iterations to convergence. Please note that we require two PCIe transfers every time because we use two different memory banks: the first transfer copies the results from one board bank (connected to the output module of the stencil accelerator) to the host memory; the second moves the results from the host to the other board bank (connected to the input module). Potentially, we could use a single one and avoid the second transfer (reducing the accelerator performance due to memory bus congestion) or overlap it with the convergence check.

If we increase the spatial parallelism from \(sp_j\) to \(sp_y\) , we can no longer place \(tp_i\) stencils due to the higher resource usage. Still, assuming no frequency changes, we can use Equation (12) to derive a threshold for temporal parallelism \(tp_x\) to ensure the application execution time does not increase:In the first approximation, we can ignore C if our application has a predefined number of iterations to execute. Similarly, we can approximate \(T_{tp_x, sp_y}\) as \(T_{1,sp_j}\cdot sp_j/sp_y\) to facilitate the computation of \(tp_x\) . In this way, it does not depend on \(tp_x\) anymore. Alternatively, we can use our performance model and explore different values for \(tp_x\) . Consequently, the previous equation changes as follows:

To exemplify, let us now take into account the Jacobi 1D single-FPGA design. Likewise Section 5.1, we generated multiple implementations for different spatial and temporal parallelism values. In particular, we considered the maximum temporal parallelism based on the suggestion by the Stencil Counter module. However, to limit the number of syntheses to run, we did not employ the stall-enabled cluster optimization, increasing the resource demand for each SST. At this point, we chose the best design with \(sp=1\) , which features 129 SSTs, as our reference design to compute the thresholds for other sp values according to equation (14). Although such an Equation does not depend on It, we selected a sufficiently large number of iterations ( \(It=1290\) ) to invoke each Jacobi 1D accelerator multiple times. Figure 14 shows how the overall execution times change according to spatial and temporal parallelism. It also reports the thresholds (vertical lines) that accurately approximate the necessary temporal parallelism to surpass the reference design’s performance. Please note that such thresholds would be different if we also considered C. Finally, it is worth noticing that the best design in terms of execution time is not \(sp=16\) , \(tp=50\) , which reaches 479.779 GFLOPS, but \(sp=8, tp=70\) (323.883 GFLOPS), proving that top results in GFLOPS do not imply top performance in a real-case scenario.

In summary, the purpose of this analysis was to highlight the trade-off between spatial and temporal parallelism, especially when considering real scenarios, which many studies tend to ignore. For this reason, we provided a formula to identify the minimum temporal parallelism threshold a stencil design has to guarantee to avoid performance degradation. Please note that the quality of the stencil implementation is a critical aspect that influences temporal parallelism. Indeed, if the design requires too many resources when the spatial parallelism grows, the resulting temporal parallelism cannot reach such a threshold. For instance, without the proposed optimizations (Section 3.2), temporal parallelism for Jacobi 1D with \(sp=4, 8\) , and 16 would have been significantly lower.

As mentioned in Section 5.2.1, comparing designs for ISLs is not trivial, even if we restrict the analysis to FPGA-based solutions only, because the aspects to consider are multiple. However, one of the main limitations to achieving such a goal is the lack of details about each literature solution (e.g., resource usage or network bandwidth), preventing the comparison of additional design features. This issue derives from the fact that most of the literature studies tend to assess the quality of their approach using GFLOPS and (sometimes) GFLOPS/W. For this reason, we expanded the range of evaluations and normalized such metrics by off-chip memory/network bandwidth, as they play a crucial role in ISL accelerators. Of course, other comparisons are viable, such as considering the semiconductor technology of each FPGA, which we reported in Tables 5 and 6. However, we are unaware of any truly effective metric for that purpose. Finally, to the best of our knowledge, we were the first to propose an exploration of trade-offs between temporal and spatial parallelism in a real-case scenario. Thus, we believe this article offers fair and adequate methods to compare Senju’s designs against the State-of-the-Art about FPGA- and SST-based ISLs.

We now discuss where top-performing designs (from Section 5.2.1) position on the Roofline Model [59]. To this end, we built Roofline Models of our target device (i.e., PAC D5005) for single- and multi-FPGA scenarios. In the former case, we used the aggregated bandwidth of the employed off-chip memory banks for the memory ceiling and the peak (32-bit) GFLOPS for the performance ceiling, computed as the number of 18 \(\times\) 19 multipliers (two per DSP on our FPGA) operating at the target frequency. Then, we plot our top-performing design results onto the Roofline Model utilizing the OI (FLOPs per transferred bytes from/to the memory banks) as the x-coordinate and the GFLOPS as the y-coordinate. In the latter case, although various approaches to building Roofline Models for a single FPGA exist in the literature [48, 49], we are unaware of any similar study for multi-FPGA systems. One way could be to consider multiple FPGAs as one single and large FPGA and aggregate their peak performance and bandwidth. However, since the network acts as the bottleneck in our multi-FPGA scenario, we propose a novel specialized formulation of the Roofline Model where we substitute the memory-bound area with a network-bound one. Please note that such a formulation is tailored to our computing scenario; thus, it may not apply to other multi-FPGA cases that do not tightly depend on the network bandwidth as we do. Specifically, we used the aggregated bandwidth of the network links for the network ceiling and the peak GFLOPS of two FPGAs for the performance ceiling. Then, to plot our multi-FPGA results, we redefined the OI as FLOPs per the aggregated amount of bytes passing through network links. For instance, in our two-FPGA scenario (chain topology), we employ one network link (100 Gbps for the network ceiling) and transfer through it \(2^{28}\cdot 4\) bytes. On the other hand, in a ring topology, we would employ two network links (200 Gbps for the network ceiling) and transfer twice the bytes ( \(2^{28}\cdot 4\) bytes per link). Consequently, the results for both topologies would be equivalent.

Figures 15(a) and 15(b) report the Roofline Models for single- and multi-FPGA scenarios, respectively. The solid black lines in each chart depict the Roofline Models of the setups we employed for the experiments, that is, 200 MHz and two memory banks for the single-FPGA case, and 200 MHz and a single network link at 100 Gbps (12.5 GB/s) for the multi-FPGA one. The colored/dashed lines indicate other setups that we will discuss later. In both scenarios, our designs are in the memory-/network-bound area of the Roofline Models, which is reasonable due to the nature of stencil computations, and reach the respective ceilings. Specifically, we are closer to the memory ceiling than the network one, thanks to better utilization of the former bandwidth (99% vs. 94%).

Table 2 reports the OI of each benchmark in a single SST execution. Thanks to its memory sub-system, an SST can buffer and reuse data through on-chip FIFOs; thus, a given SST needs to read the input data from the off-chip memory (or the previous SST) only once, limiting the impact of off-chip data transfer on the OI. Besides, temporal parallelism is particularly effective at increasing OI because it enables passing data from an SST to the next one without accessing the off-chip memory. However, even if we enough resources to make OI higher than the ridge-point (i.e., the point where the memory and performance ceilings meet), we cannot achieve higher performance than the peak one. One way to increase the performance further would be to run our designs at a higher frequency, such as 300 MHz, to fully saturate the off-chip memory bandwidth, as indicated by the solid red line in Figure 15(a). Unfortunately, we cannot reach such a frequency due to our shell limitations. Alternatively, we could exploit the remaining two off-chip memory banks and adopt a different approach, similar to the one by SASA [54] with six HBM banks; still, such a design choice would notably limit temporal parallelism, which is probably more important than spatial parallelism in real-case scenarios (see Section 5.3). Besides, as proved in Section 5.2.1, our approach performs better than SASA’s under different metrics employing one-third of the memory banks. Finally, one way to boost performance in the multi-FPGA scenario is leveraging the latest/future FPGA-based boards [1, 31], which reach 400 Gbps communication through 112 Gbps PAM4 transceivers and QSFP-DD (dashed colored line in Figure 15(b)). Of course, this approach is viable if the target board supplies enough off-chip memory bandwidth; otherwise, the memory would become the bottleneck.

Finally, we examine the limitations of Senju. First, Senju currently supports the design flow and performance modeling for Intel PAC D5005 only. Thus, an extension to other boards would require updating/changing the model and the system-level integration of our accelerators within a new shell. Speaking of the shell, its memory controllers limit the maximum frequency we can reach (around 250 MHz if we synthesize the shell only), preventing us from fully exploiting the capability of our board. Nonetheless, we achieved state-of-the-art results even with a lower frequency, which was in line with other literature studies [6, 42, 58]. Finally, Senju currently does not support ISLs containing spatial dependencies between grid point updates (e.g., Seidel 2D [3, 44]) or requiring multiple input/output streams (e.g., FDTD [3]). The former case is a historical limitation of SST-based designs. Indeed, previous studies either did not consider this kind of stencil [42] or implemented non-pipelined architectures [6, 40], significantly lowering the final GFLOPS. The latter case can represent a good opportunity to exploit more off-chip memory banks. Nonetheless, we believe that Senju already offers valuable features and performance that can be extended in future studies.