An Adaptive Sparse Matrix-Vector Multiplication Strategy Selection and Optimization Method

Abstract

The present invention discloses an adaptive sparse matrix-vector multiplication strategy selection and optimization method, which is applicable to the GPU architecture. The method includes: dividing the matrix to be processed by rows, counting the number of non-zero elements in each matrix sub-block, and if the difference multiple of the number of non-zero elements in each matrix sub-block is higher than the first preset threshold, then using the adaptive CSR-vector algorithm for processing; counting the average number of non-zero elements per row of the matrix to be processed, and if the average number of non-zero elements per row of the matrix is lower than the second preset threshold, then using the improved CSR-stream algorithm for solution; counting the number of non-zero elements of the matrix to be processed, and if the number of its non-zero elements is greater than the third preset threshold, then using the hola algorithm for solution; if none of the above conditions are met, then using the CSR-vector algorithm for solution. The present invention realizes the adaptive and efficient SpMV solution for different application problems.

Description

Technical Field

[0001] The present invention relates to the field of high-performance computing technology, and particularly relates to an adaptive sparse matrix-vector multiplication strategy selection and optimization method applicable to a GPU architecture. Background Art

[0002] SpMV (Sparse matrix–vector multiplication) is a matrix-vector multiplication operation in the form of y = αAx + βy (where A is a sparse matrix, x and y are dense vectors, and α and β are scalars), and it has wide applications in fields such as scientific computing, economic construction, signal processing, and document retrieval. Usually, the matrix involved in the calculation is sparse, such as the matrix generated by discretizing a physical process. In recent years, SpMV has been classified as one of the seven numerical methods considered to be at least very important for science and engineering in the next decade.

[0003] For different storage formats of the sparse matrix A, there will be different SpMV algorithms. The most widely used sparse matrix storage format is the CSR format. Currently, various optimizations have been carried out on the SpMV algorithm for the CSR storage format based on heterogeneous multi-core computing systems in the world, and its performance has been significantly improved. For example, the CSR-Vector algorithm, the CSR-Stream algorithm, and the hola algorithm. The algorithm ideas of the three algorithms are as follows:

[0004] 1) CSR-Vector algorithm: The wavefront is divided into multiple thread vectors (vectors) containing an equal number of threads, where one vector contains multiple threads. Finally, one vector calculates the multiplication of one row of the matrix and the x vector, and one thread in the vector is responsible for calculating part of the data in one row. After one vector finishes the calculation, the results of each thread are reduced to the 0th thread, and the 0th thread in the vector writes the result back to the y array in the memory. The disadvantage of CSR-Vector is also obvious. For the case where the number of non-zero elements in a matrix row is small, when performing the steps of loading data from memory and writing back data, the memory access between vectors is discontinuous. However, for the case where the number of non-zero elements in a matrix row is large, this disadvantage can be effectively alleviated.

[0005] 2) CSR-Stream algorithm: It divides according to the number of rows of the matrix, and each block is responsible for calculating the same number of rows. In the CSR-Stream method, each thread within a block continuously accesses the value array and the column index array of non-zero elements Col_Index from memory, and saves their multiplication results in the LDS. Finally, several threads perform reduction from the LDS and each thread reduces one row. In the CSR-Stream method, the accesses of each thread to the Value array storing non-zero elements, the column index array of non-zero elements Col_Index, and the y vector are all continuous. However, since this method needs to pre-allocate storage space in the LDS for the product of the elements in the Value array and the x vector, the CSR-Stream method is not applicable to matrices with too many non-zero elements per row.

[0006] 3) hola algorithm: It divides tasks according to the number of non-zero elements, so that each block calculates the same number of non-zero elements. Since in the CSR format, the non-zero values and column indices of the matrix are arranged in row order, it is also more suitable to adopt "task division according to non-zero elements", especially for matrices with a relatively large average number of non-zero elements per row, where the load / store overhead of its y array and row offset index Row_Index is relatively small. Using "the number of non-zero elements" for load balancing is a very effective task division strategy. The disadvantage of this method is that it requires a preprocessing process on the GPU side, which may not bring obvious performance improvement for small-scale matrices, and may even lead to performance degradation.

[0007] For the above algorithms, although they improve the computing performance under specific matrix conditions, their generality for solving problems in different fields (such as structural mechanics, computational fluid dynamics, computer graphics, etc.) is relatively poor. Summary of the Invention

[0008] The present invention provides an adaptive sparse matrix-vector multiplication strategy selection and optimization method to solve the technical problem that although the existing SpMV algorithms improve the computing performance under specific matrix conditions, their generality for solving problems in different fields is relatively poor.

[0009] To solve the above technical problems, the present invention provides the following technical solutions:

[0010] On the one hand, the present invention provides an adaptive sparse matrix-vector multiplication strategy selection and optimization method applicable to the GPU architecture. The adaptive sparse matrix-vector multiplication strategy selection and optimization method includes:

[0011] Block the matrix to be processed by rows to obtain multiple matrix sub-blocks;

[0012] Count the number of non-zero elements in each matrix sub-block. If the difference multiple of the number of non-zero elements in each matrix sub-block is higher than the first preset threshold, then use the adaptive CSR-Vector algorithm to process each matrix sub-block;

[0013] Count the average number of non-zero elements per row of the matrix to be processed. If the average number of non-zero elements per row of the matrix to be processed is lower than the second preset threshold, then use the improved CSR-Stream algorithm to solve it;

[0014] Count the number of non-zero elements of the matrix to be processed. If the number of non-zero elements of the matrix to be processed is greater than the third preset threshold, then use the hola algorithm to solve it;

[0015] If none of the above conditions are met, then use the CSR-Vector algorithm to solve it.

[0016] Further, using the adaptive CSR-Vector algorithm to process each matrix sub-block includes:

[0017] Taking wavefront as the unit of resource allocation, for each matrix sub-block, allocate different numbers of wavefronts for calculation according to the number of non-zero elements;

[0018] Inside the wavefront, set different VECTOR_SIZE for different matrix sub-blocks.

[0019] Further, the calculation method of VECTOR_SIZE is:

[0020]

[0021] where nnz_block i is the number of non-zero elements of the i-th matrix sub-block, row i is the number of rows of the i-th matrix sub-block, represents rounding down, and Pow(x) represents the smallest integer power of 2 that is greater than or equal to x.

[0022] Further, using the improved CSR-Stream algorithm to solve it includes:

[0023] Initialization, including: Given positive integers rows_per_block and block_lds_size, representing the number of rows processed by a block and the size of the shared memory owned by the block respectively; Set the number of blocks and the number of threads within the block, where the number of blocks is m is the number of rows of the matrix to be processed; the number of threads within the block is a fixed value; Set the parameter VECTOR_SIZE according to the average number of non-zero elements per row:

[0024]

[0025] Inside the kernel function: Calculate the range of rows processed by the current block, and determine the starting row and the ending row;

[0026] According to the starting row and the ending row of the block, obtain the starting index and the ending index of the non-zero element array and the column index array that the current block needs to access;

[0027] According to the starting index and the ending index of the non-zero element array and the column index array that the block needs to access, calculate the number of non-zero elements that the current block needs to process;

[0028] If the number of non-zero elements that the block needs to process is not greater than block_lds_size, then call the CSR_Stream algorithm to calculate the matrix rows in the current block;

[0029] If the number of non-zero elements that the block needs to process is greater than block_lds_size, then use the CSR-Vector idea for calculation. Specifically: First, according to the number of rows n_block_rows that the current block actually needs to calculate, the number of vectors in the block, and the VECTOR_SIZE parameter, divide the n_block_rows rows among the vectors, and each vector executes the calculation process of the standard CSR-Vector algorithm.

[0030] The beneficial effects brought by the technical solution provided by the present invention at least include:

[0031] 1) Use a wider type of sparse matrix, analyze the characteristics of the sparse matrix from multiple indicators, and classify the sparse matrix in more detail.

[0032] 2) Adaptive solution strategy selection. Without changing the efficiency of the original solution algorithm, different SpMV solution algorithms are selected for sparse matrices with different characteristics, solving the problem of high-efficiency solution of a single SpMV algorithm that cannot adapt to different sparse matrix characteristics under a general GPU, and realizing the high-efficiency SpMV solution of sparse matrices with different matrix characteristics generated for different solution problems under a general GPU architecture.

[0033] 3) Adaptive setting of solution parameters. According to the matrix characteristics, different solution parameters (such as VECTOR_SIZE) are set for the solution algorithm (such as the CSR-Vector algorithm), more deeply combining the matrix characteristics with the characteristics of the GPU architecture, and achieving higher SpMV algorithm solution efficiency. BRIEF DESCRIPTION OF THE DRAWINGS

[0034] To more clearly illustrate the technical solutions in the embodiments of the present invention, the following will briefly introduce the accompanying drawings required for the description of the embodiments. Obviously, the accompanying drawings in the following description are only some embodiments of the present invention. For those of ordinary skill in the art, without creative efforts, other accompanying drawings can be obtained based on these drawings.

[0035] Figure 1 It is a schematic execution flowchart of an adaptive sparse matrix-vector multiplication strategy selection and optimization method provided by an embodiment of the present invention;

[0036] Figure 2 It is a schematic diagram of a sparse matrix with a block property provided by an embodiment of the present invention;

[0037] Figure 3 It is a schematic diagram of a "short row matrix" provided by an embodiment of the present invention;

[0038] Figure 4 It is a schematic diagram of a "large matrix" provided by an embodiment of the present invention. Specific Embodiments

[0039] To make the objectives, technical solutions, and advantages of the present invention clearer, the following will further describe the embodiments of the present invention in detail with reference to the accompanying drawings.

[0040] First Embodiment

[0041] Since various SpMV algorithms have their own characteristics and their effects in different scenarios are also different, in order to maximize the efficiency of various SpMV algorithms and achieve better application effects, this embodiment provides an adaptive sparse matrix-vector multiplication strategy selection and optimization method applicable to the GPU architecture. This method is an adaptive strategy selection (Adaptive) algorithm that selects the most appropriate SpMV algorithm based on the characteristics of the matrix. It involves three existing SpMV algorithms: CSR-Vector algorithm, CSR-Stream algorithm, and hola algorithm. This Adaptive algorithm realizes the analysis of matrix characteristics, and the time consumption of this part can be completely ignored. It selects the most appropriate calculation strategy according to characteristics such as the average number of non-zero elements per row of the matrix, the total number of non-zero elements of the matrix, and the size of the matrix. Its process is as Figure 1 shown and includes:

[0042] S1, Matrix Blocking

[0043] The matrix A to be processed is blocked by rows to obtain multiple matrix sub-blocks; among them, each matrix sub-block contains the same or approximately the same number of rows.

[0044] S2, Determine Whether to Calculate by Block

[0045] Count the number of non-zero elements in each matrix sub-block. If the number of non-zero elements in each matrix sub-block varies greatly, such as Figure 2 shown (for example, the number of non-zero elements in the first 1 / 2 rows of the matrix is several times that of the last 1 / 2 rows), then: 1) Use the adaptive CSR-Vector algorithm (Adaptive CSR-Vector algorithm) for processing, that is, use the wavefront as the unit of resource allocation, and allocate different numbers of wavefronts for calculation according to the corresponding number of non-zero elements; 2) Inside the wavefront, set different VECTOR_SIZE for different matrix sub-blocks. Specifically, when dividing the vector in the CSR-Vector method, set different numbers of threads for the vector according to the average number of non-zero elements per row of different matrix blocks. Among them, Figure 2 the black dots in it represent the non-zero elements in the matrix. As Figure 2 shown, there are a large number of long rows (many non-zero elements per row) in the first half of the matrix. If this matrix is evenly divided into two blocks by row, then the number of non-zero elements contained in the first 1 / 2 rows of the matrix is several times that of the last 1 / 2 rows. At this time, the matrix is recognized as a sparse matrix with block properties.

[0046] Specifically, in this embodiment, the implementation process of the above steps is as follows:

[0047] Assume that the size of the sparse matrix A is m╳n, then the CSR format of A (all indexes start from 0) is:

[0048] Value[]: Used to store the values of the non-zero elements in the matrix, with a length equal to the number of non-zero elements;

[0049] Col_Index[]: The i-th element records the column number of the Value[i] element, with a length equal to the number of non-zero elements;

[0050] Row_Index[]: The i-th element records the number of non-zero elements contained in the first i-1 rows, with a length equal to the number of rows of the matrix.

[0051] Assume that the matrix A is divided into two blocks bp0 = Row_Index[m / 4], bp1 = Row_Index[m / 2], bp2 = Row_Index[3m / 4], bp3 = Row_Index[m], then the number of non-zero elements in the two matrix sub-blocks is nnz_block_0 = bp1 - 0, nnz_block_1 = bp3 - bp1;

[0052] Judge the difference in the number of non-zero elements between two matrix sub-blocks. Assume that if the difference in the number of non-zero elements in each block is more than k times, it is defined as a large gap. Then, when nnz_block_0 / nnz_block_1 ≥ k or nnz_block_1 / nnz_block_0 ≥ k, the Adaptive CSR-Vector algorithm is adopted. For two matrix sub-blocks, the CSR-Vector algorithm with different configuration parameters is used for solution according to the number of non-zero elements. And according to the average number of non-zero elements per row of different matrix sub-blocks, different VECTOR_SIZE parameters are set for the vector. Specifically: where, nnz_block i is the number of non-zero elements of the i-th matrix block, and row i is the number of rows of the i-th matrix block. represents rounding down, and Pow(x) represents the smallest integer power of 2 that is greater than or equal to x.

[0053] S3. Determine whether it is a "short-row matrix"

[0054] Count the average number of non-zero elements per row of the matrix to be processed. Taking the average number of non-zero elements per row of the matrix as an index, a matrix with an average number of non-zero elements lower than a certain value is defined as a "short-row matrix". As Figure 3 shown, if matrix A is determined to be a "short-row matrix", the improved CSR-Stream algorithm is used for solution.

[0055] where, Figure 3 the black dots in it represent non-zero elements in the matrix. As Figure 3 shown, the number of non-zero elements in each row of the matrix is very small. If the average number of non-zero elements is less than the set threshold (set by the user according to the specific solution problem, and setting it to 5 according to experience can achieve higher efficiency), it is determined as a "short-row matrix".

[0056] Specifically, in this embodiment, the implementation process of the above steps is as follows:

[0057] Calculate the average number of non-zero elements per row avg_nnz_per_row = bp3 / m. Assume that a matrix with an average number of non-zero elements less than short_row_num is defined as a "short-row matrix". Then, when avg_nnz_per_row < short_row_num, the "improved CSR-Stream algorithm" is adopted. The specific steps are as shown in S3.1~S3.4, where S3.2~S3.6 are the calculation processes of the kernel function.

[0058] S3.1, Initialization: (1) Given positive integers rows_per_block and block_lds_size, representing the number of rows processed by a block and the size of the shared memory owned by the block respectively. (2) Set the number of blocks and the number of threads within a block, where the number of blocks is and the number of threads within a block is a fixed value. (3) In the same way as in S2 above, set the parameter VECTOR_SIZE according to the average number of non-zero elements per row.

[0059] S3.2, Inside the kernel function: Obtain the id of the current block, denoted as block_id. Then calculate the range of rows processed by the current block, with the starting row id being block_row_begin and the ending row id being block_row_end (excluding this ending row). Among them, block_row_begin = block_id * rows_per_block, block_row_end = min(block_row_begin + rows_per_block, m).

[0060] S3.3, According to the starting and ending rows of the block, obtain the starting index and ending index of the non-zero element array Col_Index and the column index array Value that the block needs to access from the Row_Index array on the device side (the accessed non-zero elements do not include this ending index). Among them, the starting index block_row_idx_begin = Row_Index[block_row_begin], and the ending index block_row_idx_end = Row_Index[block_row_end].

[0061] S3.4, Calculate the number of non-zero elements that the block needs to process, block_non_zeros = block_row_idx_end - lock_row_idx_begin, based on block_row_idx_begin and block_row_idx_end.

[0062] S3.5, If block_non_zeros ≤ block_lds_size, call the CSR_Stream method to calculate the matrix rows in this block.

[0063] S3.6, if block_non_zeros > block_lds_size, the CSR-Vector idea is adopted for calculation. Specifically: First, according to the number of rows n_block_rows actually required for calculation in this block (n_block_rows = block_row_end - block_row_begin), the number of vectors in the block, and the VECTOR_SIZE parameter, divide the n_block_rows rows among the vectors, and perform the calculation process of the standard CSR-Vector algorithm within each vector.

[0064] S4. Determine whether it is a "large matrix"

[0065] Count the number of non-zero elements in the matrix to be processed. Taking the number of non-zero elements in the matrix as an index, define the matrix with the number of non-zero elements greater than a certain value (large_matrix_num) as a "large matrix". For example, Figure 4 As shown, if matrix A is determined to be a "large matrix", the hola algorithm is used for solution. That is, when bp3 > large_matrix_num, the hola algorithm is used for solution.

[0066] Among them, Figure 4 The black dots in represent the non-zero elements in the matrix. For example, Figure 4 In the shown matrix, the total number of non-zero elements is relatively large. If the total number of non-zero elements is greater than the set threshold (set by the user according to the specific solution problem, and set to 2 according to experience 23 can achieve higher efficiency), it is determined to be a "large matrix".

[0067] S5. If none of the above conditions are met, the non-block CSR-Vector algorithm is used for solution.

[0068] In summary, the method of this embodiment first divides the matrix into multiple blocks, then judges the differences between the block matrices. If the differences are too large, the Adaptive CSR-Vector algorithm is used. Otherwise, it is judged whether it is a "short row matrix". If so, the CSR-Stream algorithm is used. Otherwise, it is judged whether it is a "large matrix". If so, the hola algorithm is used. Otherwise, the CSR-Vector algorithm is used. Thus, it combines the matrix characteristics of the sparse matrix and the characteristics of the GPU architecture, adaptively selects the solution algorithm based on the sparse matrix characteristics, and adaptively sets the corresponding solution parameters for the solution algorithm, thereby realizing the adaptive and efficient SpMV solution for different application problems.

[0069] It should also be noted that, in this text, relational terms such as first and second are only used to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any actual relationship or order between these entities or operations. The term "comprising", "including" or any other variant thereof is intended to cover non-exclusive inclusion, so that a process, method, article or terminal device comprising a series of elements not only includes those elements, but also includes other elements not explicitly listed, or further includes elements inherent to such process, method, article or terminal device. Without further limitation, an element defined by the statement "comprising an..." does not exclude the presence of additional identical elements in the process, method, article or terminal device comprising the said element.

[0070] Finally, it should be noted that the above description is the preferred embodiment of the present invention. It should be pointed out that although the preferred embodiments of the present invention have been described, for those skilled in the art of this technology, once the basic creative concept of the present invention is known, several improvements and refinements can be made without departing from the principle of the present invention, and these improvements and refinements should also be regarded as the protection scope of the present invention. Therefore, the appended claims are intended to be interpreted as including the preferred embodiments and all changes and modifications falling within the scope of the embodiments of the present invention.

Claims

1. An adaptive sparse matrix-vector multiplication strategy selection and optimization method, applicable to the GPU architecture, characterized in that The adaptive sparse matrix-vector multiplication strategy selection and optimization method includes: The matrix to be processed is block-divided by rows to obtain multiple matrix sub-blocks; Count the number of non-zero elements in each matrix sub-block. If the difference multiple of the number of non-zero elements in each matrix sub-block is higher than the first preset threshold, the adaptive CSR-Vector algorithm is used to process each matrix sub-block; Count the average number of non-zero elements per row of the matrix to be processed. If the average number of non-zero elements per row of the matrix to be processed is lower than the second preset threshold, the improved CSR-Stream algorithm is used for solution; Count the number of non-zero elements of the matrix to be processed. If the number of non-zero elements of the matrix to be processed is greater than the third preset threshold, the hola algorithm is used for solution; If none of the above conditions are met, the CSR-Vector algorithm is used for solution; The use of the improved CSR-Stream algorithm for solution includes: Initialization, including: given positive integers rows_per_block and block_lds_size, representing the number of rows processed by a block and the size of the shared memory owned by the block respectively; setting the number of blocks and the number of threads within a block, where the number of blocks is m is the number of rows of the matrix to be processed; the number of threads within a block is a fixed value; setting the parameter VECTOR_SIZE according to the average number of non-zero elements per row: where nnz_block i is the number of non - zero elements of the i - th matrix sub - block, row i is the number of rows of the i - th matrix sub - block, represents rounding down, and Pow(x) represents the smallest integer power of 2 that is greater than or equal to x; Inside the kernel function: Calculate the range of rows processed by the current block to determine the start row and the end row; According to the start row and the end row of the block, obtain the start index and the end index of the non-zero element array and the column index array to be accessed by the current block; According to the start index and the end index of the non-zero element array and the column index array to be accessed by the block, calculate the number of non-zero elements to be processed by the current block; If the number of non-zero elements to be processed by the block is not greater than block_lds_size, call the CSR_Stream algorithm to calculate the matrix rows in the current block; If the number of non-zero elements to be processed by the block is greater than block_lds_size, the CSR-Vector idea is used for calculation. Specifically: First, according to the number of rows n_block_rows actually to be calculated by the current block, the number of vectors in the block, and the VECTOR_SIZE parameter, divide the n_block_rows rows among the vectors, and the standard CSR-Vector algorithm calculation process is executed within each vector.

2. The adaptive sparse matrix-vector multiplication strategy selection and optimization method according to claim 1, characterized in that The use of the adaptive CSR-Vector algorithm to process each matrix sub-block includes: Taking the wavefront as the unit of resource allocation, for each matrix sub-block, different numbers of wavefronts are allocated for calculation according to the number of non-zero elements; Inside the wavefront, different VECTOR_SIZEs are set for different matrix sub-blocks.

Images