A matrix data architecture and its method for accelerating smvm based on compressed sparse column

Abstract

The application relates to a matrix data architecture, comprising a matrix head and a matrix block; the matrix head comprises a matrix row, a column number, a non-zero element number, a matrix block number and a vectorization length S VVec ; the matrix block comprises a matrix block element; the length of an x vector pixel geometry block corresponding to each of the matrix blocks is S ImgB ; the structure of the matrix block is at least two, which are a fixed-length dense VxG block and a fixed-length sparse VxG block. A method for accelerating sparse matrix vector multiplication (SPMV) based on compressed sparse columns is used to convert the matrix format of an input sparse matrix A into the above-mentioned matrix data architecture, which is a sparse matrix storage format; the SPMV calculation is performed on a local temporary y vector reordering and full vectorization. The data architecture has the characteristics of small space occupation. In combination with the method for accelerating SPMV based on compressed sparse columns, hardware vectorization calculation instructions and memory bandwidth can be efficiently utilized, so that the SPMV can run at a faster speed, and the purpose of acceleration is achieved.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of optimization acceleration, and more particularly, to a matrix data architecture and a method for accelerating SPMV based on compressed sparse columns. BACKGROUND

[0002] Sparse matrix vector multiplication (SPMV) is a core routine used in many applications: scientific and engineering computations, economic modeling, graph analysis based on social network data, etc. In these applications, the SPMV operation is usually performed with the same matrix at a high frequency and dominates the performance. Therefore, accelerating SPMV on various specific hardware, especially on a single computer node, is an important task. In a single computer node, the central processing unit (CPU) cores share the main memory to exchange data and are the main general-purpose computing devices.

[0003] The performance of SPMV on a single computer node is limited by the memory bandwidth and the instruction latency in the computation. The memory bandwidth determines the shortest time to transfer matrix data between the processor and the memory. When using the full computer node, the peak performance of SPMV is determined by the memory bandwidth. The instruction latency usually comes from the floating-point instructions and the instructions to load and store vector elements. In addition, for SPMV, the latency also comes from the instructions to indicate the non-zero elements. When using a small number of CPU threads, the time consumed by SPMV is mainly determined by the instruction latency, in which the latency generated by the instructions to indicate the non-zero elements is negligible. There are two types of methods to accelerate SPMV generated by integral equations, such as iterative CT image reconstruction. The first type of method utilizes the sparse pattern of a specific matrix and obtains ideal acceleration in specific CT image reconstruction algorithms. The second type of method ignores the sparse pattern of a specific matrix, converts the matrix into a general storage format, and uses the commonly used matrix storage format to optimize the SPMV operation, such as the compressed sparse row (CSR), the compressed sparse column (CSC), the block-wise format, and the hybrid format (SPC5, LAV), etc.

[0004] In order to improve the efficiency of using wide vector instructions, some acceleration SPMV calculation methods use local dense matrix block to represent sparse matrix, but the effect is closely related to the application. The reason is that after using local dense matrix block, a large number of useless 0 elements are added to the matrix, causing useless memory operation and calculation operation. The acceleration SPMV calculation method uses a global vector and matrix rearrangement scheme, and before calculating SPMV once, the vector is reordered once, and after the calculation is completed, the vector is reordered again. The difficulty lies in that most of the parameters required for reordering are determined on a case-by-case basis, and the calculation time for determining the parameters is even tens to thousands of times of the time for performing SPMV once, which greatly affects the practical application of the acceleration SPMV calculation method. SUMMARY

[0005] The application provides a matrix data architecture and an acceleration SPMV method based on compressed sparse columns, which can efficiently utilize hardware vectorization calculation instructions and memory bandwidth, thereby accelerating SPMV.

[0006] To solve the above technical problems, the technical scheme adopted by the application is as follows: a matrix data architecture, comprising a matrix header and a matrix block; the matrix header comprises the number of matrix rows and columns, the number of non-zero elements, the number of matrix blocks, and a vectorization length S VVec ; the length of the x vector pixel geometric block corresponding to each matrix block is S ImgB ; the structure of the matrix block has at least two types, which are fixed-length dense VxG blocks and fixed-length sparse VxG blocks. VxG is the abbreviation of Vectorized eXecution Group.

[0007] The data architecture is obtained by processing the format of the sparse matrix, and has the characteristics of occupying less space, and the use of two different matrix block structures can achieve different degrees of compression of storage. When the matrix data architecture is used to accelerate SPMV calculation, the matrix block structure parameters can be quickly determined.

[0008] Preferably, the fixed-length dense VxG block and the fixed-length sparse VxG block each comprise a matrix block header and a VxG array; the matrix block header comprises the row and column range description of the sub-matrix, the number of VxGs, the number of columns S VxG contained in each VxG, and the number of matrix block elements contained in each column S E ; the VxG array comprises a matrix block element array, an index array of x components corresponding to each VxG, and a range of corresponding y vectors; the non-zero element dense sub-matrix corresponding to the i th VxG in the matrix block is stored in the S VxG *(i-1)*SE to the S VxG *(i)*S E -1 element, and the matrix block elements are stored in a row-block-first manner.

[0009] Preferably, the non-zero elements in the same column are stored in multiple matrix block elements, and the dense matrix block element is a vector with a fixed length of S VVec , and the matrix elements in the same dense matrix block element CSCVE have the same column index and consecutive row indices.

[0010] Preferably, the VxG array of the fixed-length sparse VxG block further comprises a non-zero element indication bitmap vector of the matrix block elements in each VxG. When the fixed-length sparse VxG block is used, an additional expansion operation is required to fill the dense matrix block element.

[0011] Preferably, the non-zero elements in the same column are stored in multiple matrix block elements, and the sparse matrix block element is a vector with a row index range length not greater than S VVec , and the matrix elements in the same sparse matrix block element have the same column index.

[0012] Preferably, the VxG connects multiple column matrix block elements sharing the same matrix row index, and ensures that the matrix row indices to be accessed in adjacent matrix block elements are the same or adjacent; the VxG corresponds to a dense sub-matrix, and the column indices of different columns on the VxG are not necessarily consecutive.

[0013] A method for accelerating SPMV based on compressed sparse columns, comprising the following steps:

[0014] Step one: converting the input sparse matrix A into a matrix format of the above-mentioned matrix data architecture, which is a sparse matrix storage format;

[0015] Step two: local temporary y vector reordering and fully vectorized SPMV calculation.

[0016] Through the local reordering strategy and the above-mentioned matrix data architecture, the method for accelerating SPMV can efficiently utilize hardware vectorization calculation instructions and memory bandwidth, so that SPMV can run at a faster speed and achieve the purpose of acceleration.

[0017] Preferably, the specific process of step two is as follows:

[0018] S1: input the sparse matrix A (size m*n) converted in step one, the vector x with a length of n, the vector y with a length of m, and the vectorization length S VVec ; input the local row vector rearrangement mapping set {ι k} and inverse mapping set {ι k -1};

[0019] S2: performing local temporary vector reordering and fully vectorized SPMV computation, such that z = local temporary vector reordering (ι k ,y) and ~z = fully vectorized SPMV (A k ,x, z ,S VVec );

[0020] S3:

[0021] S4: taking y~ in S3 as y vector in S2, performing fully vectorized SPMV computation on each sub-matrix block A k of matrix A. k Looping steps S2 and S3 until all sub-matrix blocks A K of matrix A are finished.

[0022] The reason that this method can optimize SPMV code and get the advantage of low-level SIMD instruction without using assembly code is that the fully vectorized SPMV algorithm can be easily optimized by the compiler and generate the required SIMD instruction code.

[0023] Preferably, in the S2, the specific process is as follows:

[0024] S2.1: performing vectorized operation q = q + E * x on the dense matrix E stored in VxG in the matrix block A K , where the initial value of q is a 0 vector with length m, and the process includes:

[0025] S2.1.1 performing vectorized operation on the i-th matrix block element M i in E to obtain W = a * M i + q i , where W, M i and q i are all vectors with length S VVec , a is a scalar which is the x component value of the i mod S VxG th element index in the x component index array of E, q i is the q vector element array corresponding to M E , whose index range is from (the starting value of the y vector corresponding to E * (i-1) * S VVec / S VxG ) to (the starting value of the y vector corresponding to E * i * S VVec / S VxG -1), and q i = W is written into the above index range of the q vector.

[0026] S2.1.2: cyclically execute step S2.1.1, let q of step S2.1.1 in the last cycle as q of step S2.1.1 in the next cycle, q in the next cycle = a*M in the last cycle i +q i ;

[0027] S2.1.3: complete vectorization operation of all matrix block elements M i in E, get ~z = q = q + E*x

[0028] S2.2: cyclically execute steps S2.1, S2.1.1, S2.1.2 and S2.1.3, let q of step S2.1.3 in the last cycle as q of step S2.1 in the next cycle, q in the next cycle = q in the last cycle = q + E*x;

[0029] S2.3: complete vectorization operation of all dense matrix E stored in VxG in A K , get ~z = fully vectorized SPMV(A k , x, z, S VVec ) = z + A k *x = z + q.

[0030] Each floating point operation directly reads a matrix element vector of a given length S vvec , then multiplies an x vector element and aligns and sums and stores to a subscript-continuous and length S vvec sub-vector on a y vector. Therefore, the calculation process of SPMV is a fully vectorized process, and there is no additional access to non-continuous vector element operation. And S vvec is a variable parameter, which can adapt to different width and hardware architecture vectorization instructions without changing the calculation process.

[0031] Preferably, the local row vector rearrangement mapping set changes the y vector according to the i, j coordinate axis ordering layout on the projection data region to the ordering layout parallel and perpendicular to the track of a pixel point in the projection data region; and the non-zero elements of each column of the matrix are block arranged and stored according to the y vector rearrangement mapping.

[0032] Compared with the prior art, the data architecture has the characteristics of occupying less space, and can speed up the operation process. In combination with the method for accelerating SPMV based on compressed sparse column, hardware vectorization calculation instruction and memory bandwidth can be efficiently utilized, so that SPMV can run at a faster speed, and the purpose of acceleration is achieved. BRIEF DESCRIPTION OF DRAWINGS

[0033] Figure 1 is a flow chart of a method for accelerating SPMV based on compressed sparse columns of the present application;

[0034] Figure 2 is the distribution of different pixel points and their trajectories in the projection region in the CT image reconstruction process;

[0035] Figure 3 is Figure 1 the CSCVE memory layout of the matrix column corresponding to different pixels;

[0036] Figure 4 is the SIMD efficiency of the vector y and the different memory layouts of the matrix non-zero elements;

[0037] Figure 5 is the full process of the single-threaded vectorized sparse matrix vector multiplication algorithm based on compressed sparse columns;

[0038] Figure 6 is the statistics of memory usage, floating point operations per second, and read-only memory bandwidth utilization when all CPU cores on 10 software platforms perform SPMV. DETAILED DESCRIPTION

[0039] The drawings are only used for illustrative description and cannot be understood as a limitation of the present patent; in order to better illustrate the embodiments, some components in the drawings may be omitted, enlarged or reduced, and do not represent the actual product size; for those skilled in the art, it is understandable that some well-known structures and their descriptions in the drawings may be omitted. The positional relationship described in the drawings is only used for illustrative description and cannot be understood as a limitation of the present patent.

[0040] The technical solutions of the present application will be further described in detail below through specific embodiments and in conjunction with the drawings:

[0041] Embodiment 1

[0042] As Figure 1 shown is an embodiment of a matrix data architecture, including a matrix header and a matrix block; the matrix header includes the number of matrix rows, columns, non-zero elements, the number of matrix blocks, and the vectorization length S VVec ; the edge length of the x vector pixel geometry block corresponding to each matrix block is S ImgB ; the structure of the matrix block has at least two types, which are fixed-length dense VxG blocks and fixed-length sparse VxG blocks. VxG is the abbreviation of Vectorized eXecution Group, which is the basic unit participating in the execution of SPMV in the matrix data architecture, and its row index is continuous, so that the storage and access of the data used for indexing can be greatly compressed.

[0043] The fixed-length dense VxG block and the fixed-length sparse VxG block each include a matrix block header and a VxG array; the matrix block header includes a row and column range description of a sub-matrix, a number of VxGs, a number of columns S contained by each VxG VxG , a number of matrix block elements contained by each column S E ; the VxG array includes a matrix block element array, an index array of x components corresponding to each VxG, and a range of y vectors; a non-zero element dense sub-matrix corresponding to the i-th VxG in the matrix block is stored in the S VxG *(i-1)*S E -th to S VxG *(i)*S E -1-th elements of the matrix block element array, and the matrix block elements are stored in a row block priority manner.

[0044] Specifically, non-zero elements on the same column are stored in multiple matrix block elements, the dense matrix block element is a vector with a fixed length of S VVec , and the matrix elements stored in the same dense matrix block element have the same column index and consecutive row indexes. The VxG array of the fixed-length sparse VxG block further includes a non-zero element indication bitmap vector of the matrix block elements in each VxG. When the fixed-length sparse VxG block is used, an additional expansion operation is required to fill the matrix block elements in a dense manner. Non-zero elements on the same column are stored in multiple matrix block elements, the sparse matrix block element is a vector with a row index range length not greater than S VVec , and the matrix elements in the same sparse matrix block element have the same column index.

[0045] In the embodiment, the VxG connects multiple column matrix block elements sharing the same matrix row index, and ensures that the matrix row indexes to be accessed in adjacent matrix block elements are the same or adjacent; the VxG corresponds to a dense sub-matrix, and the column indexes of different columns on the VxG are not necessarily consecutive.

[0046] The matrix format with the matrix data architecture form of the embodiment is referred to as a CSCV format matrix, and the matrix block element is referred to as a CSCVE. The compressed sparse column vector is referred to as a CSCV.

[0047] The working principle of the embodiment is that the data architecture is obtained by processing the format of a sparse matrix, and the data architecture has the characteristics of occupying less space, and the use of two different matrix block structures can achieve different degrees of compression of storage. When the matrix data architecture is used to accelerate SPMV calculation, the matrix block structure parameters can be quickly determined.

[0048] The embodiment has the advantages that the data architecture has the characteristics of occupying less space, and the running speed can be accelerated in the process of program operation.

[0049] Embodiment 2

[0050] A method for accelerating SPMV based on compressed sparse column, comprising the following steps:

[0051] Step one: converting the matrix format of the input sparse matrix A into the matrix format of the matrix data architecture of Embodiment 1;

[0052] Step two: SPMV calculation with local temporary y vector reordering and full vectorization.

[0053] Through the local reordering strategy and the above-mentioned matrix data architecture, the method for accelerating SPMV can efficiently utilize hardware vectorization calculation instructions and memory bandwidth, so that SPMV can run at a faster speed and achieve the purpose of acceleration.

[0054] Specifically, the specific process of step two is as follows:

[0055] S1: input the sparse matrix A (size m*n) converted in step one, the vector x with a length of n, the vector y with a length of m, and the vectorization length S VVec ; input the local row vector rearrangement mapping set {ι k} and the inverse mapping set {ι k -1};

[0056] S2: perform SPMV calculation with local temporary y vector reordering and full vectorization on the sub-matrix block A k of the matrix A, so that z = local temporary vector reordering (ι k , y) and = full vectorization SPMV (A k , x, z , S VVec ); = full vectorization SPMV (A k , x, z , S VVec ) has the following specific process:

[0057] S2.1: perform vectorization operation q=q+E*x on the dense matrix E stored in VxG in the matrix block A K , the initial value of q being a 0 vector with a length of m, the process including

[0058] S2.1.1 perform vectorization operation on the i-th matrix block element M i in E to obtain W=a*M i +q i , wherein W, M i , and qi are vectors of length S VVec , a is a scalar, and a VxG is the x-component value of the i mod S i th element of the subscripted array of x-components corresponding to E, and q E is the corresponding q-vector element array with subscript range from (the starting value of the y-vector corresponding to E * (i-1) * S VVec / S VxG ) to (the starting value of the y-vector corresponding to E * i * S VVec / S VxG -1), and q i = W is written into the above subscript range of the q-vector.

[0059] S2.1.2: Loop S2.1.1, let q of step S2.1.1 in the last loop be q of step S2.1.1 in the next loop, and q in the next loop = a * M i + q i in the last loop;

[0060] S2.1.3: Complete the vectorization of all matrix block elements M i in E, get ~z = q = q + E * x

[0061] S2.2: Loop S2.1, S2.1.1, S2.1.2, and S2.1.3, let q of step S2.1.3 in the last loop be q of step S2.1 in the next loop, and q in the next loop = q = q + E * x in the last loop;

[0062] S2.3: Complete the vectorization of all dense matrices E stored in VxG in A K , get = fully vectorized SPMV (A k , x, z, S VVec ) = z + A k * x = z = q.

[0063] S3: Let vector = local temporary vector reordering

[0064] S4: Loop S2 and S3 for each sub-matrix block A k of matrix A with as the y-vector in S2 until all sub-matrix blocks A k complete the fully vectorized SPMV computation.

[0065] This method directly reads a given length Svvec matrix element vector is then multiplied by an x vector element and aligned summed and stored into a sub vector of y vector with contiguous index and length S vvec The process of SPMV computation is a fully vectorized process without extra access non-contiguous index vector element operation. Moreover, S vvec is a variable parameter which can make the computation process adaptable to different width and hardware architecture vectorization instruction without changing the process.

[0066] The reason that this method can optimize SPMV code without using assembly code and get the advantage of low level SIMD instruction is that the fully vectorized SPMV algorithm can be easily optimized by compiler automatically and generate the required SIMD instruction code.

[0067] The beneficial effects of this embodiment: the matrix format of the matrix data architecture form of embodiment 1 combined with the method of accelerating SPMV based on compressed sparse column can efficiently utilize hardware vectorization calculation instruction and memory bandwidth, so that SPMV can run at a faster speed and achieve the purpose of acceleration.

[0068] Embodiment 3

[0069] Embodiment 3 of a method of accelerating SPMV based on compressed sparse column, based on embodiment 2, is used for image reconstruction process, and the local vector rearrangement mapping (referred to as IOBLR) of the line integral operator equation discretization matrix is to arrange the y vector according to the i, j coordinate axes on the projection data region to become the arrangement according to the trajectory parallel and perpendicular to a pixel point in the projection data region. The non-zero elements of each column of the matrix are arranged and stored according to the y vector rearrangement mapping. The arrangement of the matrix is called IOBLR priority layout. The inverse mapping of IOBLR is to rearrange the y vector into the layout parallel to the i, j coordinate axes on the projection data region.

[0070] The display of the selected IOBLR in this embodiment on the projection data is shown in Figure 3 . Figure 3 The arrow broken line of Figure 3 shows the IOBLR used for the arrangement of the y vector, and the square on the same arrow broken line shows the result of the rearrangement of the matrix elements under the IOBLR. Figure 1Memory layout of a 3-column sub-matrix before and after conversion to the matrix data architecture of Embodiment 1. Before conversion to the matrix data architecture, the non-zero elements of each column correspond to the points on the pixel projection trajectory, represented by the blue grid. The non-zero elements in the column are stored in the probe cell-first layout. After conversion, the non-zero elements recorded in the matrix data architecture are represented by the grid on the red parallel broken line. After conversion, the newly added padding zeros are represented by the yellow grid. Each matrix block element (CSCVE) corresponds to an arrow broken line segment. The non-zero values in the matrix block element (CSCVE) are stored in the IOBLR-first layout, in which the elements in the arrow direction are continuous in memory. The shape of the parallel broken line is determined by the minimum number of probe cells curve of the reference pixel, which is determined by the center point of the pixel block. Figure 1 Different colors in the figure represent different pixel points, and yellow and gray represent the trajectories common to multiple pixel points. The pixel points represent the geometric distribution of the x vector elements in the image reconstruction problem in SPMV, and the projection data represent the geometric distribution of the y vector elements in the image reconstruction problem in SPMV, and the trajectories of the same pixel point represent the geometric distribution of the non-zero elements of the same column matrix in the image reconstruction problem in SPMV. Figure 2 The blue and yellow grids in the figure represent the non-zero elements and zero elements in the system matrix, respectively; the colored grids in the arrow broken line constitute the CSCVE. The broken line shows the trajectory of the reference pixel in the IOBLR.

[0071] In this embodiment, the display of the selected IOBLR on the projection data is as shown in Figure 3 , which includes the rearrangement mapping of the row vector y and the memory layout of the matrix elements under the rearrangement of the row vector (referred to as IOBLR-first storage layout). In Figure 4 , the three different memory layouts: probe cell-first, angle-first, and IOBLR-first, are represented by arrow lines of different colors. The SIMD efficiency is defined by the area of the intersection of the arrow line and the blue grid, where the blue grid represents the non-zero elements in the matrix. When S VVec is 8, the ranges of the probe cell-first, angle-first, and IOBLR-first are 3, 2-6, and 7-8, respectively.

[0072] IOBLR provides efficient y-vector data local reordering and matrix element reordering schemes for compressed sparse column based SPMV algorithms and the ability to quickly determine the matrix block structure parameters. The effectiveness of IOBLR stems from the geometric properties of the integral operator in the image reconstruction process. This property exists in imaging matrices of different sizes. In addition, selecting local IOBLR can obtain an effective ordering scheme without post-analysis. Therefore, the selection of CSCV matrix block structure parameters does not require case-by-case analysis; as long as an imaging matrix is selected according to the imaging geometry of each hardware platform and a small number of parameter combinations SPMV are measured, the running time can be obtained by comparing the running time.

[0073] In this embodiment, the software implementation of the compressed sparse column based SPMV algorithm includes two branches: CSCV-Z and CSCV-M. The compressed sparse column based SPMV algorithm can be adjusted in combination with other optimization methods, such as the removal of 0-element filling technology. The calculation processes of the two implementations are only slightly different, that is, whether to remove the 0-element filling in CSCVE. CSCV-M uses the 0-element filling removal technology, which can further reduce the memory usage, and the vexpend operation needs to be performed during the reading process of the matrix CSCVE; while the memory bandwidth usage of CSCV-Z is higher, and no additional matrix element processing is required. The input matrix block of CSCV-M uses the fixed-length sparse VxG block format, while the input matrix block of CSCV-Z uses the fixed-length dense VxG block format.

[0074] In this embodiment, the single-threaded calculation process (including local temporary reordering and fully vectorized SPMV calculation in the calculation process) of the two implementations of CSCV-Z and CSCV-M is as shown in Figure 5

[0075] In order to demonstrate the performance advantage of the CSCV SPMV calculation process, this embodiment selects two x86 architecture hardware platforms and 10 current mainstream SPMV software implementations for comparative testing. The case used for testing is to reconstruct an image of 1024*1024 pixels according to the parallel beam geometry CT on all CPU cores on the platform. The test process only performs forward projection calculation, and the matrix size is 700800*1048576, and the number of non-zero elements is about 1.3 billion. The test hardware platforms include: Intel Xeon Gold 6130 (referred to as SKL) and AMD EPYC 7452 (referred to as Zen2). The vector instruction width of the two platforms is different. The SKL platform uses 512-bit width vectorization instructions, while the Zen2 platform uses 256-bit width vectorization instructions. The number of CPU cores of the two platforms is different, SKL is 32, and Zen2 is 64.

[0076] ​The SPMV software implementations used in the tests of this embodiment include two software implementations of CSCV, CSCV-Z and CSCV-M, and SPC5, CVR, Merge, CSR5, VHCC, ESB, Intel MKL CSC (referred to as MKL-CSC) and CSR (referred to as MKL-CSR). Two data precisions, single precision and double precision, are included in the tests. Because the clinical CT imaging of hospitals is all single precision, this kind of test case is more important.

[0077] The tests of this embodiment mainly investigate three aspects: memory usage of matrix storage, floating point operation per second and read-only memory bandwidth utilization. The specific test data is shown in Table 2. Figure 6 According to the test data, the memory usage M Rit of CSCV-M is reduced to 53% to 68% of that of MKL-CSR, the floating point operation per second Performance (GFLOPS) can reach 3.51 times (i.e. the calculation speed is 3.51 times) of that of MKL-CSR, and the read-only memory bandwidth utilization R EM is greater than 85% and the highest is 98%. The top two of the floating point operation per second are the software implementations of CSCV. Compared with the third of the experiment, the floating point operation per second is 1.6 to 2 times that of the latter. For short word length data (single precision data), the speedup is more obvious.

[0078] In order to show that CSCV SPMV has the advantages of adapting to wide SIMD instructions and obtaining low-level SIMD instructions without using assembly code, we compare the floating point operation per second of CSCV-M and SPC5. In the tests of this embodiment, CSCV-M does not use assembly code, and the floating point operation per second is obviously 1.6 times that of SPC5 which uses bottom assembly code, and can be used across hardware platforms. The SPMV calculation kernel of SPC5 is written in vectorized calculation process using Intel AVX-512 assembly code and manually optimized, and can only run on Intel's platform. CSCV-M is only written in high-level language C++, and the vectorization process is realized by the automatic optimization of the compiler. The CSCV SPMV calculation process has a full vectorization calculation kernel, which is easy to be optimized by the compiler, so it still has good performance, and can be used across hardware platforms.

[0079] In order to show that the calculation process of the method of this embodiment has the advantage of quickly determining the block structure parameters, this embodiment introduces four different size CT imaging matrices to perform the same performance test as before. The parameters of the matrices are shown in Table 1. The software implementation of CSCV uses the same parameters as before (see Table 2), including the vector length S vvec , the size of the block S ImgB and the length of VxG SVxG The statistics of the floating point operation per second of various software in the four matrix tests are shown in Table 3. The results show that the CSCV software implementation has obvious computational advantage, and the maximum value of floating point operation per second is superior to that of the opponent.

[0080] Table 1 Parameters of the four CT imaging matrices to be tested

[0081]

[0082] Table 2 Block architecture parameters adopted by different software implementations in running on different platforms

[0083]

[0084]

[0085] Table 3

[0086]

[0087] Obviously, the above embodiments of the present application are only examples for clearly illustrating the present application, and are not intended to limit the implementation manners of the present application. Based on the above description, other different forms of changes or variations can be made by those skilled in the art. Here, it is not necessary and also impossible to enumerate all the implementation manners. Any modification, equivalent replacement and improvement, etc. made within the spirit and principle of the present application shall be included in the protection scope of the claims of the present application.

Claims

1. A method for accelerating SPMV based on compressed sparse columns, characterized in that, Includes the following steps: Step 1: Convert the input sparse matrix A from its matrix format into a matrix data architecture, which is a sparse matrix storage format. The matrix data architecture includes a matrix header and matrix blocks. The matrix header includes the number of matrix rows and columns, the number of non-zero elements, the number of matrix blocks, and the vectorization length S. VVec The side length of the x-vector pixel geometric cube corresponding to each of the matrix blocks is S. ImgB The matrix block has at least two structures: a fixed-length dense VxG block and a fixed-length sparse VxG block. Both the fixed-length dense VxG block and the fixed-length sparse VxG block include a matrix block header and a VxG array; the matrix block header includes a description of the row and column ranges of the submatrix, the number of VxGs, and the number of columns S contained in each VxG. VxG The number of matrix block elements in each column, S E The VxG array includes a matrix block element array, an array of indices corresponding to the x-components of each VxG, and a range of the corresponding y-vectors; the non-zero dense submatrix corresponding to the i-th VxG in the matrix block is stored in the S-th element array of the matrix block. VxG ×(i-1)×S E Up to the Sth VxG ×(i)×S E -1 elements, and store matrix block elements in row-matrix order; Step 2: Reorder the local temporary y vectors and compute the fully vectorized SPMV.

2. The method for accelerating SPMV based on compressed sparse columns according to claim 1, characterized in that, In the matrix data architecture, non-zero elements in the same column are stored in multiple matrix block elements, and dense matrix block elements have a fixed length. A vector whose matrix elements stored in the same dense matrix block have the same column index and consecutive row index.

3. The method for accelerating SPMV based on compressed sparse columns according to claim 1, characterized in that, The VxG array of the fixed-length sparse VxG block also includes a bitmap vector indicating the non-zero elements of each VxG matrix block element.

4. The method for accelerating SPMV based on compressed sparse columns according to claim 3, characterized in that, In the matrix data architecture, non-zero elements in the same column are stored in multiple matrix block elements, and sparse matrix block elements are a range of row indices with a length not exceeding a certain value. A vector in which matrix elements in the same sparse matrix block have the same column index.

5. The method for accelerating SPMV based on compressed sparse columns according to claim 1, characterized in that, The VxG connection shares access to multiple matrix block elements with the same matrix row index, and ensures that the matrix row index to be accessed in adjacent matrix block elements is the same or adjacent; the VxG corresponds to a dense submatrix and the column indices of different columns in the VxG are not necessarily consecutive.

6. A method for accelerating SPMV based on compressed sparse columns according to any one of claims 1-5, characterized in that, The specific process for step two is as follows: S1: Input the sparse matrix A after format conversion in step one, a vector x of length n, and a vector y of length m, with a vectorized length of S. VVec Input local row vector rearrangement mapping set { } and the set of inverse maps { }; S2: Perform local temporary y-vector reordering and fully vectorized SPMV computation, so that... =Local temporary vector reordering( ,y) and =Fully vectorized SPMV(A k ,x, ,S VVec ); S3: Let vector y = =Local temporary vector reordering( , ); S4: Based on S3 As the y-vector in S2, for each submatrix block A of matrix A k Repeat steps S2 and S3 until all submatrix blocks A are reached. k Complete the fully vectorized SPMV computation.

7. The method for accelerating SPMV based on compressed sparse columns according to claim 6, characterized in that, In S2, the specific process is as follows: S2.1: For matrix block A K The dense matrix E stored in VxG is vectorized into q = q + E × x, where the initial value of q is a zero vector of length m. The process includes: S2.1.1 For the i-th matrix block element M in E i Vectorization operation yields W = a × M i +q i Among them, W and M i and q i All are of length S VVec The vector, where 'a' is a scalar and is the i-th modulo S in the array of subscripts corresponding to the x-components of E. VxG The x-component value of each element index, q i It is M E The corresponding array of elements of the q vector has an index range from (E to the starting value of the y vector × (i-1) × S) VVec / S VxG (E corresponds to the initial value of the y vector × i × S) VVec / S VxG -1), and q i =W is written into the above index range of the q vector; S2.1.2: Repeat step S2.1.1, setting q from the previous iteration as q in the next iteration, where q in the next iteration equals a × M from the previous iteration. i +q i ; S2.1.3: Complete all matrix block elements M in E. i Vectorization operations yield =q=q+E×x S2.2: Repeat steps S2.1, S2.1.1, S2.1.2 and S2.1.3, and let q in step S2.1.3 in the previous loop be used as q in step S2.1 in the next loop. The q in the next loop = the q in the previous loop = q + E × x; S2.3: Complete A K Vectorization operations on all dense matrices E stored in VxG yield the following results: =Fully vectorized SPMV(A k ,x, ,S VVec )=z+A k ×x = z + q.

8. The method for accelerating SPMV based on compressed sparse columns according to claim 1, characterized in that, The local row vector rearrangement mapping set transforms the y-vector sorting layout according to the i and j coordinate axes on the projected data region into a sorting layout according to the trajectory of a certain pixel in the projected data region, which is parallel and perpendicular to the pixel. The non-zero elements of each column of the matrix are arranged and stored in blocks according to the y-vector rearrangement mapping.

Images