Data Processing Method and Device for SIMD-Facing Parallel LU Direct Solution

Abstract

The present application provides a data processing method and apparatus for parallel LU direct solution oriented to SIMD. The method includes: obtaining data to be processed, including 2m pattern matrices, each pattern matrix having the same non-zero element distribution pattern; performing a real-number processing on each pattern matrix to obtain a corresponding real-number matrix for each pattern matrix; performing a permutation on the real-number matrix corresponding to each pattern matrix to obtain a corresponding transformation matrix and a corresponding permutation matrix, where the transformation matrix is a matrix obtained after the real-number matrix is transformed, and the permutation matrix is a matrix reflecting the transformation process of the real-number matrix; performing a vectorization process based on the transformation matrices and permutation matrices corresponding to the 2m pattern matrices to obtain a one-dimensional array for storage; using the SIMD extension component to read the vectorized one-dimensional array, performing parallel LU decomposition and direct solution, and outputting the direct solution result. Thus, parallel linear solution is performed through the SIMD extension component to achieve accelerated solution.

Description

Technical Field

[0001] This application relates to the technical field of data processing. Specifically, it relates to a data processing method and device for parallel LU direct solution oriented to SIMD. Background Art

[0002] In the field of numerical analysis (such as the field of numerical analysis of electromagnetic fields), there is often a need to solve matrices. Traditional solution schemes include direct solution methods, iterative solution methods, etc.

[0003] The direct solution method usually requires LU decomposition (or called LU factorization). LU decomposition (or called LU factorization) means that given an n×n matrix A, if two matrices L and U can be found such that: A = LU, where L is a unit lower triangular matrix, that is, a lower triangular matrix with all 1s on the diagonal, and U is an upper triangular matrix, then A is said to be LU decomposable. LU decomposition is a method of decomposing a matrix into the product of a lower triangular matrix (L) and an upper triangular matrix (U). This technique is very useful for solving linear equations because it can simplify the calculation process and is usually more efficient than directly using Gaussian elimination.

[0004] Of course, the upper triangular matrix U can also be obtained by Gaussian elimination, and at the same time, record the transformation matrix corresponding to each row operation. The product of these transformation matrices is the inverse matrix of L. Therefore, by appropriately adjusting this process, L can be directly obtained.

[0005] Doolittle's method: This method assumes that L is a unit lower triangular matrix (diagonal elements are 1) and directly calculates the elements of L and U according to the formula.

[0006] Crout's method: Similar to the Doolittle method, but it assumes that U is an upper triangular matrix with diagonal elements of 1.

[0007] Cholesky decomposition: When A is a symmetric positive definite matrix, it can be further simplified to Cholesky decomposition, where A = LLT.

[0008] Once the LU decomposition of matrix A is obtained, solving the linear equation system in the form of Ax = b becomes simple. Ax = b can be replaced by LUx = b. Then, let Ux = y, so first solve Ly = b, and then solve Ux = y. Both of these two steps involve forward substitution and backward substitution of triangular systems, and both are relatively fast operations.

[0009] However, not all matrices can be LU decomposed. For example, if matrix A is singular or near - singular (i.e., its determinant is close to zero), then there may not exist such L and U such that A = LU.

[0010] Modern processor architectures (the mainstream ones include the x86 architecture and the arm architecture, etc. ARM NEON is a SIMD instruction set under the ARM platform, and the x86 platform has instruction sets such as MMX, SSE, AVX, etc. Some other architectures also have some SIMD instruction sets) can support the execution of multiple instructions simultaneously through their SIMD extension components (vectorization components). SIMD, that is, single instruction multiple data, means that a single instruction stream can operate on multiple data streams, thereby improving the operation speed. The SIMD extension component refers to the dedicated hardware unit in the processor for executing SIMD operations. These extension components usually include wide registers, dedicated SIMD instruction sets, and corresponding execution units. The following are several common SIMD extension components and their characteristics.

[0011] 1. Intel's SIMD extension: Name: SSE (Streaming SIMD Extensions); Versions: SSE, SSE2, SSE3, SSSE3, SSE4.1, SSE4.2; Register width: 128 bits; Data types: single - precision floating - point numbers (32 bits), double - precision floating - point numbers (64 bits), integers (8 bits, 16 bits, 32 bits, 64 bits); Uses: Widely used in floating - point operations, multimedia processing, and scientific computing; Registers: XMM0 to XMM15 (possibly more on newer processors).

[0012] 2. AMD's SIMD extension: Name: 3DNow!; Versions: 3DNow!, 3DNow! Professional, 3DNow! +, Enhanced 3DNow!; Register width: 64 bits; Data types: single - precision floating - point numbers (32 bits); Uses: Mainly used in floating - point operations, especially in games and multimedia applications; Registers: MM0 to MM7.

[0013] 3. ARM's SIMD extension: Name: NEON registers; Widths: 64 bits and 128 bits; Data types: single - precision floating - point numbers (32 bits), double - precision floating - point numbers (64 bits), integers (8 bits, 16 bits, 32 bits, 64 bits); Uses: Widely used in multimedia processing, signal processing, and high - performance computing; Registers: Q0 to Q31 (128 bits), D0 to D31 (64 bits).

[0014] Taking NEON as an example, the NEON technology is an advanced SIMD (Single Instruction Multiple Data) architecture for Arm Cortex-A series processors. It can accelerate multimedia and signal processing algorithms, such as video encoders / decoders, 2D / 3D graphics, games, audio and speech processing, image processing, telephony, and sound. NEON instructions perform "packed SIMD" processing, and registers are considered vectors of elements of the same data type. The supported data types are: signed / unsigned 8-bit, 16-bit, 32-bit, and 64-bit.

[0015] There are a total of 32 NEON vector registers in ARMv8, with a length of 128 bits. They can be used for processing scalar operands or vector operands. These registers can store multiple operands of the same data type, such as 2 double-type operands of 64 bits or 4 float-type operands of 32 bits. During the execution of NEON instructions, multiple operands on the same vector register are processed simultaneously, thus achieving parallel processing of different data. To accelerate the program using vector registers, the data in memory needs to be read into them in advance. Therefore, the continuity of the accessed data in memory has an important impact on the performance of vector registers, and data alignment can significantly improve the performance of vector registers.

[0016] In some fields represented by electromagnetic computing, there is usually a need to solve systems of equations with the same distribution positions of non-zero elements. For example, in the field of electromagnetic computing, when the model and mesh discretization remain unchanged, the distribution positions of non-zero elements in the systems of equations generated at different frequencies are the same. Since there are a large number of complex matrices (dense complex matrices, sparse complex matrices, etc.) that need to be linearly solved in the field of electromagnetic numerical analysis, but the current technology requires calling a solution function each time for solving, and it is impossible to perform parallel solution of 2 or more linear equations.

[0017] Based on this, the inventors of the present application designed a parallel LU direct solution scheme for SIMD to achieve accelerated solution. Summary of the Invention

[0018] The purpose of the embodiments of the present application is to provide a data processing method and device for parallel LU direct solution for SIMD, so as to perform parallel linear solution calculations through SIMD extension components and achieve accelerated solution.

[0019] To achieve the above purpose, the embodiments of the present application are implemented as follows:

[0020] In a first aspect, the embodiments of the present application provide a data processing method for parallel LU direct solution for SIMD, including: obtaining data to be processed, where the data to be processed includes 2 m mode matrices, m ∈ Z+ , each pattern matrix has the same non - zero element distribution pattern; if the pattern matrix is a complex matrix, perform a real - number conversion process on each pattern matrix to obtain the real - number matrix corresponding to each pattern matrix; perform a permutation on the real - number matrix corresponding to each pattern matrix to obtain the corresponding transformation matrix and the corresponding permutation matrix, where the transformation matrix is the matrix obtained after the real - number matrix is transformed, and the permutation matrix is the matrix reflecting the transformation process of the real - number matrix; based on the transformation matrices and permutation matrices corresponding to 2 m pattern matrices, perform a vectorization process to obtain a vectorized one - dimensional array for storage; use the SIMD extension component to read the vectorized one - dimensional array, perform parallel LU decomposition and direct solution, and output the direct solution result.

[0021] Combined with the first aspect, in the first possible implementation manner of the first aspect, performing a real - number conversion process on each pattern matrix to obtain the real - number matrix corresponding to each pattern matrix includes: for each pattern matrix, perform the real - number conversion process using the following formula:

[0022] Kx = b,

[0023]

[0024] K is the current pattern matrix, x is the unknown vector, b is the right - hand - side vector, K r and K i are respectively the real part and the imaginary part of the pattern matrix K, b r and b i are respectively the real part and the imaginary part of the right - hand - side vector b, x r and x i are respectively the real part and the imaginary part of the unknown vector x, K r , - K i , - K i , - K r correspond to four sub - matrices in the form of real - number matrices, denoted as the first sub - matrix, the second sub - matrix, the third sub - matrix, and the fourth sub - matrix; distribute the first sub - matrix, the second sub - matrix, the third sub - matrix, and the fourth sub - matrix in the form from left to right and from top to bottom to form the real - number matrix corresponding to the pattern matrix.

[0025] Combined with the first aspect, in the second possible implementation manner of the first aspect, performing a permutation on the real - number matrix corresponding to each pattern matrix to obtain the corresponding transformation matrix and the corresponding permutation matrix includes: for the real - number matrix A j corresponding to each pattern matrix, j ∈ [1, 2 m : create a permutation matrix P j , and multiply both sides of A j x = b j by P j, there are:

[0026] P j A j x j = P j b j ,

[0027] Among them, the permutation matrix P j is initially the identity matrix, Aj is the real matrix corresponding to the pattern matrix, and x j is the vector of unknowns, and b j is the right-hand vector; for each column in the real matrix A j , determine the element with the largest absolute value from the column where the pivot position is located in the order of pivots, swap the row containing this element with the row where the pivot position is located, and update the permutation matrix P j to record these row exchanges, obtaining the transformation matrix A j corresponding to the real matrix A j′ and the corresponding permutation matrix P j , there is A j′ x j = P j b j , where the pivot represents the element on the main diagonal of the matrix, and the pivot position represents the position of the element on the main diagonal of the matrix.

[0028] Combined with the second possible implementation manner of the first aspect, in the third possible implementation manner of the first aspect, based on the transformation matrix and the permutation matrix corresponding to 2 m pattern matrices, perform vectorization processing to obtain a vectorized one-dimensional array for storage, including: combining the elements at the same position in all the transformation matrices A j′ corresponding to the pattern matrices into a group of vectors corresponding to this position, obtaining a vectorized one-dimensional array A for storage; combining the elements at the same position in all the right-hand vectors P j b j corresponding to the pattern matrices into a group of vectors corresponding to this position, obtaining a vectorized one-dimensional array B for storage.

[0029] Combined with the third possible implementation manner of the first aspect, in the fourth possible implementation manner of the first aspect, use the SIMD extension component to read the vectorized one-dimensional array and perform parallel LU decomposition and direct solution, including: using the SIMD extension component to read the vectorized one-dimensional array A and perform parallel LU decomposition to obtain the lower triangular matrix L j , upper triangular matrix U j ; combining all the lower triangular matrices L jThe elements at the same position are combined into a set of vectors corresponding to that position, and the vectorized one-dimensional array L is obtained for storage. Also, the upper triangular matrix U corresponding to all pattern matrices j The elements at the same position are combined into a set of vectors corresponding to that position, and the vectorized one-dimensional array U is obtained for storage; the SIMD extension component reads the vectorized one-dimensional arrays L, U, and B, and performs parallel direct solution calculations.

[0030] Combined with the fourth possible implementation of the first aspect, in the fifth possible implementation of the first aspect, the SIMD extension component reads the vectorized one-dimensional array A, performs parallel LU decomposition, and obtains the lower triangular matrix L corresponding to each pattern matrix j and the upper triangular matrix U j , including: initialization is of double-precision type and the value is set to 0.0; the SIMD extension component reads each set of vectors in the one-dimensional array A one by one, each set of vectors contains 2 m elements, and corresponding parallel multiply-accumulate calculations are performed based on the 2 m elements in each set of vectors to update values, and finally the updated lower triangular matrix L corresponding to each pattern matrix is obtained j and the upper triangular matrix U j .

[0031] Combined with the fourth possible implementation of the first aspect, in the sixth possible implementation of the first aspect, the SIMD extension component reads the vectorized one-dimensional arrays L, U, and V, and performs parallel direct solution calculations, including: based on the lower triangular matrix L corresponding to each pattern matrix j and the upper triangular matrix U j , there are:

[0032]

[0033] Initialize each intermediate vector y j , j ∈ [1, 2 m ; the SIMD extension component reads each set of vectors in the one-dimensional array L and each set of vectors in the one-dimensional array B one by one, each set of vectors contains 2 m elements, and corresponding parallel multiply-accumulate calculations are performed based on the 2 m elements in each set of vectors to solve each intermediate vector y j ; and the SIMD extension component reads each set of vectors in the one-dimensional array U one by one, each set of vectors contains 2 m elements, and based on the 2 m elements in each set of vectors and each intermediate vector yj Perform corresponding parallel multiply-accumulate calculations to solve each unknown vector x j .

[0034] In a second aspect, an embodiment of the present application provides a data processing device for parallel LU direct solution oriented to SIMD, including: a data acquisition unit for acquiring data to be processed, where the data to be processed includes 2 m mode matrices, m ∈ Z + , and each mode matrix has the same non-zero element distribution pattern; a real-number conversion unit for performing real-number conversion processing on each mode matrix when the mode matrix is a complex matrix to obtain a real matrix corresponding to each mode matrix; a matrix permutation unit for permuting the real matrix corresponding to each mode matrix to obtain a corresponding transformation matrix and a corresponding permutation matrix, where the transformation matrix is a matrix obtained after the real matrix is transformed, and the permutation matrix is a matrix reflecting the transformation process of the real matrix; a vectorization unit for performing vectorization processing based on the transformation matrices and permutation matrices corresponding to 2 m mode matrices to obtain a vectorized one-dimensional array for storage; a parallel solution unit for using the SIMD extension component to read the vectorized one-dimensional array, perform parallel LU decomposition and direct solution, and output a direct solution result.

[0035] In a third aspect, an embodiment of the present application provides a storage medium, which is provided in an electronic device and includes a stored program. When the program runs, it controls the electronic device where the storage medium is located to execute the data processing method for parallel LU direct solution oriented to SIMD according to any one of the first aspect or possible implementation manners of the first aspect.

[0036] In a fourth aspect, an embodiment of the present application provides an electronic device, including a memory and a processor. The memory is used to store information including program instructions, and the processor is used to control the execution of the program instructions. When the program instructions are loaded and executed by the processor, the steps of the data processing method for parallel LU direct solution oriented to SIMD according to any one of the first aspect or possible implementation manners of the first aspect are implemented.

[0037] Beneficial effects:

[0038] 1. In the field of electromagnetic numerical analysis (or other numerical analysis fields, such as the field of mechanical numerical analysis), there are a large number of complex matrices (dense complex matrices, sparse complex matrices, etc.) that need to be linearly solved. In order to use the SIMD extension component to implement parallel LU direct solution and improve the solution efficiency, it is necessary to consider the conversion problem of complex matrices and the problem that LU decomposition cannot be performed. Accordingly, in this solution, by acquiring the data to be processed, including 2 m mode matrices, m ∈ Z +, each pattern matrix has the same non - zero element distribution pattern. If the pattern matrix is a complex matrix, it is made real, and a real matrix corresponding to each pattern matrix is obtained. If it is a real matrix, no real - number conversion is required, thus solving the problem of real - number conversion of complex matrices. By permuting the real matrix corresponding to each pattern matrix, the corresponding transformation matrix and the corresponding permutation matrix are obtained. The transformation matrix is the matrix obtained after the real matrix is transformed, and the permutation matrix is the matrix reflecting the transformation process of the real matrix. The permutation operation is based on determining the pivot element of the matrix. For each column in the real matrix A j , the element with the largest absolute value is determined from the column where the pivot position is located in the order of the pivot, the row containing this element with the largest absolute value is swapped with the row where the pivot position is located, and the permutation matrix P is updated j to record these row exchanges, obtaining the real matrix A j corresponding transformation matrix A j′ and the corresponding permutation matrix P j . The pivot represents the element on the main diagonal of the matrix, and the pivot position represents the position of the element on the main diagonal of the matrix. Such a permutation operation based on columns can effectively avoid the situation where LU decomposition cannot be performed or calculation errors occur because the matrix is singular or near - singular. Then, for the transformation matrices and permutation matrices corresponding to 2 m pattern matrices, vectorization is performed to obtain a vectorized one - dimensional array for storage. Then, the SIMD extension component is used to read the vectorized one - dimensional array for parallel LU decomposition and direct solution, and the direct solution result is output. Such a scheme can effectively utilize the acceleration performance of the SIMD extension component to achieve parallel linear solution calculation, greatly improving the calculation efficiency and shortening the solution time.

[0039] 2. Combine the elements at the same position of the transformation matrices A j′ corresponding to all pattern matrices into a group of vectors corresponding to this position, and obtain a vectorized one - dimensional array A for storage; combine the elements at the same position of the right - hand - side vectors P j b j corresponding to all pattern matrices into a group of vectors corresponding to this position, and obtain a vectorized one - dimensional array B for storage. In this way, the elements of all pattern matrices at the same position can be combined into a vector at this position and stored in the form of a one - dimensional array (aligned). (The SIMD extension component can store multiple operands of the same data type, such as 2 64 - bit double - type operands or 4 32 - bit float - type operands, etc.). When the SIMD extension component reads later, the jump step can be set to 2 m, so when performing computational operations, "packed SIMD" processing can be carried out using SIMD instructions, enabling multiple operands on the same vector register to be processed simultaneously, thus achieving parallel processing of different data. And using vector registers to accelerate the program requires reading the data in memory into them in advance. Data alignment can significantly improve the performance of SIMD extension components, effectively enhancing the computational speed of parallel processing, and thus improving the direct solution efficiency.

[0040] 3. During the solution process, through permutation, LU decomposition, and solution calculation, this process can be illustrated as follows: After calculating y in Ly = Pb, substitute it into Ux = y to solve for the unknown vector x. The SIMD extension component can play a role in the processes of LU decomposition, solution calculation, etc. to achieve parallel processing. For example, read the vectorized one-dimensional array A and perform parallel LU decomposition to obtain the lower triangular matrix L corresponding to each mode matrix j , upper triangular matrix U j ; Combine the elements of the lower triangular matrix L j at the same position into a group of vectors corresponding to that position, obtain the vectorized one-dimensional array L for storage, and combine the elements of the upper triangular matrix U j at the same position into a group of vectors corresponding to that position, obtain the vectorized one-dimensional array U for storage; Use the SIMD extension component to read the vectorized one-dimensional array L, one-dimensional array U, and one-dimensional array B, and perform parallel direct solution calculation. Such a method can play a great role in the solution process of numerous mode matrices (such as sparse complex matrices, dense complex matrices, etc. with the same form). It can not only use the SIMD extension component to achieve parallel processing, but also ensure the continuity of the accessed data in memory through an aligned data storage scheme, greatly improving the acceleration performance of the SIMD extension component. For example, in this solution, test cases are provided for two dense complex matrices of size 942 * 942. The comparison method is to directly solve LU for the two matrices without using the SIMD extension component and record the solution time. When directly solving LU, it takes 750.794s, while in this solution (taking NEON parallel direct solution as an example), it only takes 141.341s, and the efficiency is improved by 81.17%.

[0041] To make the above objects, features, and advantages of the present application more obvious and understandable, the following specifically gives preferred embodiments and, in conjunction with the accompanying drawings, makes a detailed description as follows. Description of the Drawings

[0042] To more clearly illustrate the technical solutions of the embodiments of the present application, the accompanying drawings required for the embodiments of the present application will be briefly introduced below. It should be understood that the following drawings only show certain embodiments of the present application and should not be regarded as limiting the scope. For those of ordinary skill in the art, other related drawings can be obtained based on these drawings without creative efforts.

[0043] Figure 1 It is a flowchart of a data processing method for directly solving parallel LU for SIMD provided by an embodiment of the present application.

[0044] Figure 2 It is the pattern matrix K 1 After being made real, it obtains the real matrix A 1 Schematic diagram of the process.

[0045] Figure 3 It is the pattern matrix K 2 After being made real, it obtains the real matrix A 2 Schematic diagram of the process.

[0046] Figure 4 It is the real matrix A 1 After performing a permutation operation, it obtains the transformed matrix A 1 ′ Schematic diagram of the process.

[0047] Figure 5 It is the real matrix A 2 After performing a permutation operation, it obtains the transformed matrix A 2 ′ Schematic diagram of the process.

[0048] Figure 6 It is for the transformed matrix A 1 ′ and the transformed matrix A 2 ′ Schematic diagram of vectorization.

[0049] Figure 7 It is a schematic diagram of combining vectors into a one-dimensional array for storage. Specific embodiments

[0050] Next, the technical solutions in the embodiments of the present application will be described in conjunction with the accompanying drawings in the embodiments of the present application.

[0051] In the field of numerical analysis (such as the field of numerical analysis of electromagnetic fields, the field of numerical analysis of mechanics, etc.), there is often a need to solve matrices. Among these matrices to be solved (including real matrices and complex matrices, such as dense complex matrices, sparse complex matrices, etc.), there are many pattern matrices. Here, the pattern matrix refers to a matrix with the same non-zero element distribution pattern. For example, in the field of electromagnetic calculation, when the model and mesh discretization remain unchanged, the positions of the non-zero elements in the equations obtained at different frequencies are the same. To solve these matrices, direct solution methods, iterative solution methods and other solution methods can be used. However, in the existing direct solution methods, a solution function needs to be called each time a solution is performed, and it is impossible to perform parallel solution of two or more linear equations.

[0052] Accordingly, this embodiment proposes a data processing method for parallel LU direct solution oriented to SIMD, which uses the SIMD extension component to achieve simultaneous horizontal solution of multiple equation systems with the same non-zero element distribution. For example, in the electromagnetic field, for a model that has been meshed discretely, when solving the electromagnetic responses at two frequencies of 100Hz and 10Hz, a system of equations needs to be solved for each frequency. The positions of the non-zero elements in these two systems of equations are exactly the same, generating two pattern matrices. Then, it is suitable to use the method proposed in the present invention to vectorize the solution process horizontally to achieve parallel solution.

[0053] Please refer to Figure 1 , Figure 1 which is the flowchart of the data processing method for parallel LU direct solution oriented to SIMD provided in this embodiment.

[0054] First, step S10 can be executed.

[0055] Step S10: Obtain the data to be processed. Among them, the data to be processed includes 2 m pattern matrices, m ∈ Z + , and each pattern matrix has the same non-zero element distribution pattern.

[0056] In this embodiment, these data to be processed, that is, the pattern matrices to be solved (usually including boundary conditions, written in the form of Kx = b), can be obtained. The number of pattern matrices processed at one time needs to meet the quantity of 2 m , m ∈ Z +, each pattern matrix has the same non - zero element distribution pattern. Currently, the latest ARM9 platform supports the simultaneous execution of 4 double - precision instructions, and it may be expanded in the future. Although the technical difficulty is high and the number of double - precision instructions that can be processed synchronously is limited (currently at the level of simultaneously executing 2 - 4 double - precision instructions), there is still great application potential. Moreover, the idea of this solution can theoretically achieve a higher level of parallel processing to greatly improve the processing efficiency. For the convenience of explanation in this embodiment, 2 pattern matrices are taken as examples for introduction, which should not be regarded as a limitation to this application.

[0057] Since the pattern matrix can be either a real - number matrix or a complex - number matrix, but the pattern matrices in the same batch (i.e., in the data to be processed) are either all real - number matrices or all complex - number matrices. For a real - number matrix, no real - number processing is required, but for a complex - number matrix, real - number processing is needed.

[0058] Therefore, if the pattern matrix is a complex - number matrix, step S20 is executed.

[0059] Step S20: If the pattern matrix is a complex - number matrix, perform real - number processing on each pattern matrix to obtain the real - number matrix corresponding to each pattern matrix.

[0060] In this embodiment, for each pattern matrix, the following formula is used for real - number processing:

[0061] Kx = b, (1)

[0062]

[0063] K is the current pattern matrix, x is the unknown vector, b is the right - hand - side vector, K r and K i are respectively the real part and the imaginary part of the pattern matrix K, b r and b i are the real part and the imaginary part of the right - hand - side vector b, x r and x i are the real part and the imaginary part of the unknown vector x, K r 、-K i 、-K i 、-K r correspond to four sub - matrices in the form of real - number matrices respectively, denoted as the first sub - matrix, the second sub - matrix, the third sub - matrix, and the fourth sub - matrix. Then, arrange the first sub - matrix, the second sub - matrix, the third sub - matrix, and the fourth sub - matrix in the form from left to right and top to bottom to form the real - number matrix corresponding to the pattern matrix. Assume that K is an n×n complex - number matrix, then the first sub - matrix, the second sub - matrix, the third sub - matrix, and the fourth sub - matrix are all n×n, and the formed real - number matrix A j is 2n×2n, where A jis the real matrix corresponding to the j-th pattern matrix, where j ∈ [1, 2 m .

[0064] For example Figure 2 and Figure 3 as shown, after real-numbering the pattern matrix K on the left side in Figure 2 , the real matrix A composed of the first sub-matrix, the second sub-matrix, the third sub-matrix, and the fourth sub-matrix on the right side is obtained 1 ; after real-numbering the pattern matrix K on the left side in 1 , the real matrix A composed of the first sub-matrix, the second sub-matrix, the third sub-matrix, and the fourth sub-matrix on the right side is obtained Figure 3 . 2 ; after real-numbering the pattern matrix K on the left side in 2 .

[0065] After obtaining the real matrix corresponding to each pattern matrix, step S30 can be further executed.

[0066] Step S30: Permute the real matrix corresponding to each pattern matrix to obtain the corresponding transformation matrix and the corresponding permutation matrix, where the transformation matrix is the matrix obtained after the real matrix is transformed, and the permutation matrix is the matrix reflecting the transformation process of the real matrix.

[0067] In this embodiment, in order to avoid the situation where the real matrix cannot be LU decomposed or calculation errors occur because the matrix is singular or nearly singular, the real matrix corresponding to each pattern matrix can be permuted to obtain the corresponding transformation matrix and the corresponding permutation matrix.

[0068] Exemplarily, for the real matrix A corresponding to each pattern matrix j :

[0069] A permutation matrix P j (initially the identity matrix) can be created. Multiply both sides of A j x = b j by P j , and we get:

[0070] P j A j x j = P j b j , (3)

[0071] where the permutation matrix P j is initially the identity matrix, A j is the real matrix corresponding to the pattern matrix, x j is the unknown vector, and b j is the right-hand vector.

[0072] The specific operations are as follows:

[0073] For each column in the real matrix A j Determine the element with the largest absolute value in the order of pivot elements from the column where the pivot element is located (excluding the elements in the row where the pivot element has been determined, that is, assuming the pivot element at (1, 1) has been determined and the pivot element at (2, 2) is currently being determined, it can only be determined from the remaining undetermined pivot element data). Swap the row containing this element with the row where the pivot element is located (if there are the same absolute values, swap with the row with the smallest row number), and update the permutation matrix P j To record these row swaps, obtaining the real matrix A j The corresponding transformation matrix A j′ And the corresponding permutation matrix P j , we have:

[0074] A j′ x j = P j b j , (4)

[0075] Among them, the pivot element represents the element on the main diagonal of the matrix, the pivot position represents the position of the element on the main diagonal of the matrix, and A j′ Is the transformation matrix after performing the permutation operation on the real matrix A j , and at this time the permutation matrix P j Is also gradually updated with the transformation operation. Here, no other characters are used to distinguish it from the initialized permutation matrix P j , which needs attention. In addition, in very few cases, there may be a situation where the pivot elements of the last row or several rows are 0 due to restricted or impossible row swaps at the end. In this case, the direct solution method cannot be used (the system can prompt that the matrix is a singular matrix and cannot perform LU decomposition), and other methods (such as the iterative method) need to be considered for solution.

[0076] As Figure 4 And Figure 5 Shown, after performing the permutation operation on the real matrix A 1 Corresponding to the pattern matrix K 1 , the corresponding transformation matrix A 1′ Is obtained, as Figure 4 ; After performing the permutation operation on the real matrix A 2 Corresponding to the pattern matrix K 2 , the corresponding transformation matrix A 2′ Is obtained, as Figure 5 .

[0077] After performing the permutation operation on the real matrix A j Corresponding to each pattern matrix, the corresponding transformation matrix A is obtainedj′ and the updated permutation matrix P j . After that, step S40 can be executed.

[0078] Step S40: Based on the transformation matrices and permutation matrices corresponding to 2 m pattern matrices, perform vectorization processing to obtain a vectorized one-dimensional array for storage.

[0079] To align the data and further improve the acceleration performance of the SIMD extension component, here, the elements of the transformation matrices A j′ corresponding to all pattern matrices at the same position can be combined into a group of vectors corresponding to this position, and a vectorized one-dimensional array A is obtained for storage.

[0080] Exemplarily, the elements of the transformation matrix A j′ corresponding to each pattern matrix at the same position can be combined into a group of vectors (including 2 m elements) corresponding to this position, and then they are arranged in the order of row and column numbers to form a one-dimensional array for storage. In this way, the continuity when the SIMD extension component reads can be ensured, and the acceleration performance can be effectively exerted.

[0081] As Figure 6 shown, the elements of the transformation matrix A 1′ and the transformation matrix A 2′ at the same position are combined into a group of vectors corresponding to this position. For example, at the position (1, 1), that is, the position of the first row and the first column, the element of the transformation matrix A 1′ at this position is -6, and the element of the transformation matrix A 2′ at this position is 3, which are combined into a group of vectors (-6, 3). The same combination is performed for other positions in the first row to determine the corresponding vectors. As Figure 7 shown (showing the form of the vectors in the first row forming an array), the two elements in the dashed box are a group of vectors (the gray background represents the value at the corresponding position of the transformation matrix A 1′ , the white background represents the value at the aligned position of the transformation matrix A 2′ and the transformation matrix A 1′ in this way, it is convenient for the subsequent SIMD extension component to read the information of two transformation matrices at one time for LU decomposition), and a vectorized one-dimensional array A is obtained in order (for example, the first row and the first column - the first row and the second column... the second row and the first column... the 2nth row and the 2nth column) for storage. It is saved using the double type, that is, a 128-bit vector can include two 64-bit double type data, which is convenient for subsequent parallel reading and calculation operations using the SIMD extension component. Here, if 4 transformation matrices are vectorized, a 128-bit vector can contain 4 32-bit float type data to adapt to the parallel processing of more matrices.

[0082] Similarly, the right-end vectors P corresponding to all pattern matrices can be j b j Combined with the elements at the same position into a group of vectors corresponding to this position, and a vectorized one-dimensional array B is obtained for storage.

[0083] Thus, the preprocessing of the data to be processed (2 m pattern matrices) is completed, and then step S50 can be executed.

[0084] Step S50: Use the SIMD extension component to read the vectorized one-dimensional array, perform parallel LU decomposition and direct solution, and output the direct solution result.

[0085] In this embodiment, the SIMD extension component can be used to read the vectorized one-dimensional array A, perform parallel LU decomposition, and obtain the lower triangular matrix L corresponding to each pattern matrix j , upper triangular matrix U j .

[0086] Due to modern processor architectures (the mainstream ones include x86 architecture and arm architecture, etc. ARM NEON is the SIMD instruction set under the ARM platform, and the x86 platform has instruction sets such as MMX, SSE, and AVX. Some other architectures also have some SIMD instruction sets), through its SIMD extension component (vectorization component), multiple instructions can be supported to be executed simultaneously. In this embodiment, NEON is taken as an example, which should not be regarded as a limitation to this application.

[0087] In order to achieve parallel processing, this embodiment makes parallel improvements on the basis of the LU direct decomposition based on the Doolittle method. First, the direct decomposition method of Doolittle will be introduced in detail.

[0088] Assume that the leading principal minors D of A i ≠0 (i = 1:n), then A = LU, denoted as:

[0089]

[0090] From matrix multiplication, we get:

[0091]

[0092] Then the formula for calculating other positions of the matrix can be expressed as follows:

[0093]

[0094] Record the L matrix and U matrix in turn through the above formula.

[0095] To achieve parallelism in this solution, it is possible to initialize as double-precision type, and set the value to 0.0; read each group of vectors in the one-dimensional array A one by one through the NEON registers, and each group of vectors contains 2 m elements. Based on the 2 m elements in each group of vectors, perform corresponding parallel multiply-accumulate calculations to update value, and finally obtain the updated lower triangular matrix L j and upper triangular matrix U j corresponding to each pattern matrix.

[0096] Specifically, the design of the decomposition algorithm (taking the parallel calculation of two matrices as an example), the code and explanations are as follows:

[0097]

[0098]

[0099]

[0100] In this way, parallel processing of LU decomposition can be achieved, and the updated lower triangular matrix L j and upper triangular matrix U j corresponding to each pattern matrix can be obtained.

[0101] After obtaining the updated lower triangular matrix L j and upper triangular matrix U j corresponding to each pattern matrix, the elements at the same position in the lower triangular matrices L j corresponding to all pattern matrices can be combined into a group of vectors corresponding to that position, and the vectorized one-dimensional array L is obtained for storage. Also, the elements at the same position in the upper triangular matrices U j corresponding to all pattern matrices can be combined into a group of vectors corresponding to that position, and the vectorized one-dimensional array U is obtained for storage. For the specific vectorized storage method, please refer to the previous text and will not be elaborated here.

[0102] After that, the SIMD extension component (such as the NEON register) can be used to read the vectorized one-dimensional arrays L, U, and B, and perform parallel direct solution calculations to output the direct solution result.

[0103] Exemplarily, based on the lower triangular matrix L j and upper triangular matrix U j corresponding to each pattern matrix, there is:

[0104]

[0105] where A j′ x j = Pj b j Denotes the substitution process, (L j U j )x j = P j b j Denotes the LU decomposition process, Denotes the subsequent solution process. During the solution process, it is necessary to first solve the intermediate vector y through L j y j = P j b j Solve the intermediate vector y j , and then substitute the intermediate vector y j into U j x j = y j , in order to solve the unknown vector x j .

[0106] Specifically, each intermediate vector y can be initialized j , j ∈ [1, 2 m , read each group of vectors in the one-dimensional array L and each group of vectors in the one-dimensional array B one by one through the NEON register. Each group of vectors contains 2 m elements, and perform corresponding parallel multiply-add calculations based on the 2 m elements in each group of vectors to solve each intermediate vector y j . For substitution, each group of vectors in the one-dimensional array U can be read one by one through the NEON register. Each group of vectors contains 2 m elements, and perform corresponding parallel multiply-add calculations with each intermediate vector y m based on the 2 j elements in each group of vectors to solve each unknown vector x j .

[0107] Specifically design the solution algorithm (taking the parallel calculation of two matrices as an example). The code and explanations are as follows:

[0108]

[0109]

[0110] Thus, each unknown vector x j , that is, the solution corresponding to each mode matrix, can be calculated. The output solution can also be further transformed and restored to the complex form. For example, since x is of dimension 2×n and the original complex matrix is of dimension n, the real part of the complex solution is equal to the element at the i-th position of x, and the imaginary part is equal to the element at the i + n-th position of x. Perform the complex restoration operation according to this rule to obtain the corresponding complex solution.

[0111] To verify this solution, the test case is two dense complex matrices of size 942×942. The comparison method is to directly solve the LU of the two matrices separately without using NEON registers, and record the solution time. The following Table 1 shows the comparison results of the two methods:

[0112] Table 1. Comparison Results of the Solution of this Embodiment and the Conventional Solution

[0113] Unit / Method NEON Parallel Direct Solver LU Direct Solver NEON Efficiency Improvement Time (s) 141.341 750.794 81.17％

[0114] It can be seen that when directly solving LU, it takes 750.794 s, while this solution (NEON parallel direct solution) only takes 141.341 s, and the efficiency is increased by 81.17%.

[0115] Based on the same inventive concept, the embodiment of the present application further provides a data processing device for parallel LU direct solution for SIMD, including:

[0116] A data acquisition unit for acquiring data to be processed, where the data to be processed includes 2 m pattern matrices, m∈Z + , and each pattern matrix has the same non-zero element distribution pattern.

[0117] A real number conversion unit for performing real number conversion processing on each pattern matrix when the pattern matrix is a complex matrix to obtain a real number matrix corresponding to each pattern matrix.

[0118] A matrix permutation unit for permuting the real number matrix corresponding to each pattern matrix to obtain a corresponding transformation matrix and a corresponding permutation matrix, where the transformation matrix is a matrix obtained after the real number matrix is transformed, and the permutation matrix is a matrix reflecting the transformation process of the real number matrix.

[0119] A vectorization unit for performing vectorization processing based on the transformation matrices and permutation matrices corresponding to 2 m pattern matrices, and storing the vectorized one-dimensional array obtained.

[0120] A parallel solution unit for using the SIMD extension component to read the vectorized one-dimensional array, perform parallel LU decomposition and direct solution, and output the direct solution result.

[0121] In this embodiment, the real number conversion unit is specifically used for: for each pattern matrix, performing real number conversion processing using the following formula:

[0122] Kx = b,

[0123]

[0124] K is the current pattern matrix, x is the unknown vector, b is the right-end vector, Kr and K i are the real and imaginary parts of the mode matrix K respectively, and b r and b i are the real and imaginary parts of the right-hand vector b respectively, and x r and x i are the real and imaginary parts of the unknown vector x respectively. K r , -K i , -K i , -K r correspond to four sub-matrices in the form of real matrices respectively, denoted as the first sub-matrix, the second sub-matrix, the third sub-matrix, and the fourth sub-matrix; the first sub-matrix, the second sub-matrix, the third sub-matrix, and the fourth sub-matrix are distributed in the form from left to right and from top to bottom to form the real matrix corresponding to the mode matrix.

[0125] In this embodiment, the matrix permutation unit is specifically used for: for each real matrix A corresponding to the mode matrix j , j ∈ [1, 2 m : create a permutation matrix P j , and multiply both sides of A j x = b j by P j , to obtain:

[0126] P j A j x j = P j b j ,

[0127] where the permutation matrix P j is initially the identity matrix, A j is the real matrix corresponding to the mode matrix, j is the unknown vector, and b j is the right-hand vector; for each column in the real matrix A j , determine the element with the largest absolute value in the column where the pivot element is located in the order of the pivot element, exchange the row containing this element with the row where the pivot element is located, and update the permutation matrix P j to record these row exchanges, to obtain the transformation matrix A j corresponding to the real matrix A j′ and the corresponding permutation matrix P j , and there is A j′ x j = P j b j , where the pivot element represents the element on the main diagonal of the matrix, and the pivot position represents the position of the element on the main diagonal of the matrix.

[0128] In this embodiment, the vectorization unit is specifically used for: based on 2m For the transformation matrix and permutation matrix corresponding to each pattern matrix, perform vectorization processing, and store the vectorized one-dimensional array, including: combining the elements at the same position of all the transformation matrices A corresponding to the pattern matrices j′ into a group of vectors corresponding to this position to obtain a vectorized one-dimensional array A for storage; combining the right-end vectors P corresponding to all the pattern matrices j b j into a group of vectors corresponding to this position at the same position to obtain a vectorized one-dimensional array B for storage.

[0129] In this embodiment, the parallel solution unit is specifically configured to: use the SIMD extension component to read the vectorized one-dimensional array A, perform parallel LU decomposition, and obtain the lower triangular matrix L corresponding to each pattern matrix j and the upper triangular matrix U j ; combine the elements of all the lower triangular matrices L corresponding to the pattern matrices j at the same position into a group of vectors corresponding to this position to obtain a vectorized one-dimensional array L for storage, and also, combine the elements of all the upper triangular matrices U corresponding to the pattern matrices j at the same position into a group of vectors corresponding to this position to obtain a vectorized one-dimensional array U for storage; use the SIMD extension component to read the vectorized one-dimensional arrays L, U, and B, and perform parallel direct solution calculations.

[0130] In this embodiment, the parallel solution unit is specifically configured to: initialize to be of double-precision type and set the value to 0.0; the SIMD extension component reads each group of vectors in the one-dimensional array A one by one, and each group of vectors contains 2 m elements, perform corresponding parallel multiply-add calculations based on the 2 m elements in each group of vectors, update the value of , and finally obtain the updated lower triangular matrix L and upper triangular matrix U corresponding to each pattern matrix j and the upper triangular matrix U j .

[0131] In this embodiment, the parallel solution unit is specifically configured to: based on the lower triangular matrix L and upper triangular matrix U corresponding to each pattern matrix j and the upper triangular matrix U j , there is:

[0132]

[0133] Initialize each intermediate vector y j , j ∈ [1, 2 m; The SIMD extension component reads each group of vectors in the one-dimensional array L and each group of vectors in the one-dimensional array B one by one. Each group of vectors contains 2 m elements, and performs corresponding parallel multiply-accumulate calculations based on the 2 m elements in each group of vectors to solve each intermediate vector y j ; And, the SIMD extension component reads each group of vectors in the one-dimensional array U one by one. Each group of vectors contains 2 m elements, and performs corresponding parallel multiply-accumulate calculations based on the 2 m elements in each group of vectors and each intermediate vector y j to solve each unknown vector x j .

[0134] This embodiment also provides a storage medium. The storage medium is disposed in an electronic device and includes a stored program. When the program runs, it controls the electronic device where the storage medium is located to execute the data processing method for parallel LU direct solution oriented to SIMD in this embodiment.

[0135] And, this embodiment also provides an electronic device, including a memory and a processor. The memory is used to store information including program instructions, and the processor is used to control the execution of the program instructions. When the program instructions are loaded and executed by the processor, the steps of the data processing method for parallel LU direct solution oriented to SIMD in this embodiment are implemented.

[0136] In summary, the embodiments of the present application provide a data processing method and device for parallel LU direct solution oriented to SIMD. In the field of electromagnetic numerical analysis (or other numerical analysis fields, such as the field of mechanical numerical analysis), there are a large number of complex matrices (dense complex matrices, sparse complex matrices, etc.) that need to be linearly solved. In order to use the SIMD extension component to implement parallel LU direct solution and improve the solution efficiency, it is necessary to consider the conversion problem of complex matrices and the problem that LU decomposition cannot be performed. Accordingly, this solution obtains the data to be processed, including 2 m pattern matrices, m ∈ Z + , and each pattern matrix has the same non-zero element distribution pattern. If the pattern matrix is a complex matrix, it is real-numbered to obtain the real matrix corresponding to each pattern matrix. If it is a real matrix, there is no need to perform real-numbering processing, so as to solve the real-numbering problem of complex matrices. By permuting each real matrix corresponding to the pattern matrix, the corresponding transformation matrix and the corresponding permutation matrix are obtained. The transformation matrix is the matrix obtained after the real matrix is transformed, and the permutation matrix is the matrix reflecting the transformation process of the real matrix. The permutation operation is based on determining the pivot element of the matrix. By performing operations on the real matrix A jFor each column, determine the element with the largest absolute value in the column where the pivot is located in the pivot order, swap the row containing this element with the row where the pivot is located, and update the permutation matrix P j to record these row exchanges, obtaining the real matrix A j The corresponding transformation matrix A j′ and the corresponding permutation matrix P j , where the pivot represents the element on the main diagonal of the matrix, and the pivot position represents the position of the element on the main diagonal of the matrix. Such a permutation operation is column-based, which can effectively avoid the situation where LU decomposition cannot be performed or calculation errors occur due to the matrix being singular or nearly singular. Then, for the transformation matrices and permutation matrices corresponding to 2 m pattern matrices, vectorization processing is performed to obtain a vectorized one-dimensional array for storage, and then the SIMD extension component is used to read the vectorized one-dimensional array for parallel LU decomposition and direct solution, and the direct solution result is output. Such a scheme can effectively utilize the acceleration performance of the SIMD extension component to achieve parallel linear solution calculation, greatly improving the calculation efficiency and shortening the solution time.

[0137] Combine the elements at the same position of the transformation matrices A corresponding to all pattern matrices j′ into a group of vectors corresponding to this position to obtain a vectorized one-dimensional array A for storage; combine the elements at the same position of the right-hand vectors P corresponding to all pattern matrices j b j into a group of vectors corresponding to this position to obtain a vectorized one-dimensional array B for storage. In this way, the elements at the same position of all pattern matrices can be combined into vectors at this position and stored in the form of a one-dimensional array (aligned). (The SIMD extension component can store multiple operands of the same data type, such as 2 64-bit double-type operands or 4 32-bit float-type operands, etc.). When the SIMD extension component reads later, the jump step can be set to 2 m , so during the calculation operation, "packed SIMD" processing can be performed using SIMD instructions, enabling multiple operands on the same vector register to be processed simultaneously, thus achieving parallel processing of different data. And using vector registers to accelerate the program requires reading the data in memory into them in advance. The alignment of the data can significantly improve the performance of the SIMD extension component, effectively improving the calculation speed of parallel processing, and thus improving the direct solution efficiency.

[0138] During the solution process, through permutation, LU decomposition, and solution calculation, this process can be illustrated as follows: After calculating y in Ly = Pb, substitute it into Ux = y to solve the unknown vector x. The SIMD extension component can play a role in processes such as LU decomposition and solution calculation to achieve parallel processing. For example, read the vectorized one-dimensional array A and perform parallel LU decomposition to obtain the lower triangular matrix L corresponding to each mode matrix j , the upper triangular matrix U j ; Combine the elements of the lower triangular matrix L j at the same position into a group of vectors corresponding to that position, obtain the vectorized one-dimensional array L for storage, and combine the elements of the upper triangular matrix U j at the same position into a group of vectors corresponding to that position, obtain the vectorized one-dimensional array U for storage; Use the SIMD extension component to read the vectorized one-dimensional array L, one-dimensional array U, and one-dimensional array B and perform parallel direct solution calculation. Such a method can play a great role in the solution process of numerous mode matrices (such as sparse complex matrices and dense complex matrices with the same form). It can not only use the SIMD extension component to achieve parallel processing, but also ensure the continuity of the accessed data in memory through an aligned data storage scheme, greatly improving the acceleration performance. For example, in this solution, test cases are provided for two dense complex matrices of size 942 * 942. The comparison method is to directly solve LU for the two matrices without using the SIMD extension component and record the solution time. When directly solving LU, it takes 750.794s, while in this solution (taking NEON parallel direct solution as an example), it only takes 141.341s, and the efficiency is improved by 81.17%.

[0139] The above are only examples of the present application and are not used to limit the protection scope of the present application. For those skilled in the art, various changes and modifications can be made to the present application. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present application shall be included in the protection scope of the present application.

Claims

1. A data processing method for SIMD-oriented parallel LU direct solution, characterized in that: include: Get the data to be processed, where the data to be processed includes 2 m pattern matrix, m∈Z + , each pattern matrix has the same distribution pattern of non-zero elements; If the pattern matrix is ​​a complex matrix, each pattern matrix is ​​converted to a real number to obtain a real number matrix corresponding to each pattern matrix; The real number matrix corresponding to each pattern matrix is ​​permuted to obtain a corresponding transformation matrix and a corresponding permutation matrix, wherein the transformation matrix is ​​a matrix obtained after the real number matrix is ​​transformed, and the permutation matrix is ​​a matrix reflecting the transformation process of the real number matrix; Based on 2 m The transformation matrix and permutation matrix corresponding to the pattern matrix are vectorized to obtain a vectorized one-dimensional array for storage; The SIMD extension component is used to read the vectorized one-dimensional array, perform parallel LU decomposition and direct solution, and output the direct solution result.

2. The data processing method for SIMD-oriented parallel LU direct solution according to claim 1, characterized in that: Each pattern matrix is ​​converted to a real number to obtain a real number matrix corresponding to each pattern matrix, including: For each pattern matrix, the following formula is used for real number processing: Kx=b, K is the current mode matrix, x is the unknown vector, b is the right end vector, K r and K i are the real and imaginary parts of the pattern matrix K, respectively, and b r and b i is the real and imaginary part of the right-hand vector b, x r and x i is the real and imaginary part of the unknown vector x, K r , -K i , -K i , -K r They correspond to four sub-matrices in the form of real matrices, which are recorded as the first sub-matrix, the second sub-matrix, the third sub-matrix, and the fourth sub-matrix; The first sub-matrix, the second sub-matrix, the third sub-matrix, and the fourth sub-matrix are distributed from left to right and from top to bottom to form a real number matrix corresponding to the pattern matrix.

3. The data processing method for SIMD-oriented parallel LU direct solution according to claim 1, characterized in that: The real number matrix corresponding to each pattern matrix is ​​permuted to obtain the corresponding transformation matrix and the corresponding permutation matrix, including: For each pattern matrix corresponding to the real matrix A j , j∈[1,2 m ]: Create a permutation matrix P j , in A j x=b j Multiply both ends by P j ,have: P j A j x j =P j b j , Among them, the permutation matrix P j Initially, it is the identity matrix, A j is the real matrix corresponding to the pattern matrix, x j is the unknown vector, b j is the right-hand vector; For a real matrix A j For each column in , determine the element with the largest absolute value from the column where the pivot position is located in the pivot order, swap the row containing the element with the row where the pivot position is located, and update the permutation matrix P j By recording these rows interchangeably, we get the real matrix A j The corresponding transformation matrix A j ′ and the corresponding permutation matrix P j , there is A j 'x j =P j b j , where the pivot represents the elements on the main diagonal of the matrix, and the pivot position represents the position of the elements on the main diagonal of the matrix.

4. The data processing method for SIMD-oriented parallel LU direct solution according to claim 3, characterized in that: Based on 2 m The transformation matrix and permutation matrix corresponding to the pattern matrix are vectorized to obtain a vectorized one-dimensional array for storage, including: The transformation matrix A corresponding to all pattern matrices j ′The elements at the same position are combined into a set of vectors corresponding to the position, and the vectorized one-dimensional array A is obtained for storage; The right-hand vector P corresponding to all pattern matrices j b j The elements at the same position are combined into a set of vectors corresponding to the position, and a vectorized one-dimensional array B is obtained for storage.

5. The data processing method for SIMD-oriented parallel LU direct solution according to claim 4, characterized in that: Use SIMD extension components to read vectorized one-dimensional arrays and perform parallel LU decomposition and direct solution, including: Use the SIMD extension component to read the vectorized one-dimensional array A and perform parallel LU decomposition to obtain the lower triangular matrix L corresponding to each pattern matrix. j , upper triangular matrix U j ; The lower triangular matrix L corresponding to all pattern matrices j The elements at the same position are combined into a set of vectors corresponding to the position, and the vectorized one-dimensional array L is stored, and the upper triangular matrix U corresponding to all pattern matrices is j The elements at the same position are combined into a set of vectors corresponding to the position, and the vectorized one-dimensional array U is obtained for storage; The SIMD extension component is used to read the vectorized one-dimensional array L, one-dimensional array U and one-dimensional array B to perform parallel direct solution calculations.

6. The data processing method for SIMD-oriented parallel LU direct solution according to claim 5, characterized in that: Use the SIMD extension component to read the vectorized one-dimensional array A and perform parallel LU decomposition to obtain the lower triangular matrix L corresponding to each pattern matrix. j , upper triangular matrix U j ,include: initialization It is of double precision type and its value is set to 0.0; The SIMD expansion unit reads each set of vectors in the one-dimensional array A one by one. Each set of vectors contains 2 m elements, based on 2 in each set of vectors m elements to perform corresponding parallel multiplication and addition calculations, and update The value of each pattern matrix is ​​finally obtained, and the updated lower triangular matrix L corresponding to each pattern matrix is ​​obtained. j , upper triangular matrix U j .

7. The data processing method for SIMD-oriented parallel LU direct solution according to claim 5, characterized in that: The SIMD extension component is used to read the vectorized one-dimensional array L, one-dimensional array U and one-dimensional array B, and perform parallel direct solution calculations, including: Based on the lower triangular matrix L corresponding to each pattern matrix j , upper triangular matrix U j ,have: Initialize each intermediate vector y j , j∈[1,2 m ]; The SIMD expansion unit reads each set of vectors in the one-dimensional array L and each set of vectors in the one-dimensional array B one by one. Each set of vectors contains 2 m elements, based on 2 in each set of vectors m elements to perform corresponding parallel multiplication and addition calculations to solve each intermediate vector y j ; And, the SIMD extension component reads each set of vectors in the one-dimensional array U one by one, each set of vectors contains 2 m elements, based on 2 in each set of vectors m elements and each intermediate vector y j Perform the corresponding parallel multiplication and addition calculations to solve each unknown vector x j .

8. A data processing device for SIMD-oriented parallel LU direct solution, characterized in that: include: A data acquisition unit is used to acquire the data to be processed, wherein the data to be processed includes 2 m pattern matrix, m∈Z + , each pattern matrix has the same distribution pattern of non-zero elements; A real number conversion unit, used for performing real number processing on each pattern matrix when the pattern matrix is ​​a complex number matrix, so as to obtain a real number matrix corresponding to each pattern matrix; A matrix permutation unit is used to permute the real number matrix corresponding to each pattern matrix to obtain a corresponding transformation matrix and a corresponding permutation matrix, wherein the transformation matrix is ​​a matrix obtained after the real number matrix is ​​transformed, and the permutation matrix is ​​a matrix reflecting the transformation process of the real number matrix; Vectorization unit for 2-based m The transformation matrix and permutation matrix corresponding to the pattern matrix are vectorized to obtain a vectorized one-dimensional array for storage; The parallel solution unit is used to read the vectorized one-dimensional array using the SIMD extension component, perform parallel LU decomposition and direct solution, and output the direct solution result.

9. A storage medium, characterized in that: The storage medium is set in an electronic device, including a stored program, wherein when the program is running, the electronic device where the storage medium is located is controlled to execute the data processing method for parallel LU direct solution for SIMD according to any one of claims 1 to 7.

10. An electronic device comprising a memory and a processor, wherein the memory is used to store information including program instructions, and the processor is used to control the execution of the program instructions, characterized in that: When the program instructions are loaded and executed by the processor, the steps of the data processing method for SIMD-oriented parallel LU direct solution described in any one of claims 1 to 7 are implemented.

Images