SIMD-oriented data processing method and device for parallel LU direct solution

Abstract

The invention provides a SIMD-oriented parallel LU direct solution data processing method and device, and the method comprises the steps: obtaining to-be-processed data which comprises 2m pattern matrixes, and each pattern matrix has the same non-zero element distribution pattern; performing real number processing on each mode matrix to obtain a real number matrix corresponding to each mode matrix; the real number matrix corresponding to each mode matrix is permuted, a corresponding transformation matrix and a corresponding permutation matrix are obtained, 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 the transformation matrixes and permutation matrixes corresponding to the 2m pattern matrixes, vectorization processing is carried out, and a one-dimensional array is obtained and stored; and reading the vectorized one-dimensional array by using an SIMD extension component, carrying out parallel LU decomposition and direct solution, and outputting a direct solution result. And parallel linear solution is carried out through the SIMD extension component, so that accelerated solution is realized.

Description

Technical Field

[0001] The present application relates to the field of data processing technology, and in particular to a data processing method and device for SIMD-oriented parallel LU direct solution. Background Art

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

[0003] Direct solutions usually require the use of LU decomposition (or LU factorization). LU decomposition (or 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 simplifies the computational process and is usually more efficient than using Gaussian elimination directly.

[0004] Of course, we can also use Gaussian elimination to obtain the upper triangular matrix U, and 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, we can directly obtain L.

[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] Krotnek-Kaplansky algorithm (Crout's method): Similar to Doolittle's method, but it assumes that U is an upper triangular matrix with all diagonal elements set to 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 the matrix A is obtained, solving a linear system of equations of the form Ax = b becomes simple. Ax = b can be replaced by LUx = b, then let Ux = y, then we can first solve Ly = b and then solve Ux = y. Both steps involve forward and backward substitution of the triangular system, both of which are relatively fast operations.

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

[0010] Modern processor architectures (mainstream ones include x86 architecture and ARM architecture, ARM NEON is the SIMD instruction set under the ARM platform, the x86 platform has MMX, SSE, AVX and other instruction sets, and some other architectures also have some SIMD instruction sets) can support the simultaneous execution of multiple instructions through their SIMD extension components (vectorization components). SIMD, which stands for single instruction multiple data, means that a single operation instruction can execute multiple data streams, thereby increasing the operation speed. SIMD extension components refer to dedicated hardware units in the processor for performing 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 extensions: name SSE (Streaming SIMD Extensions); versions: SSE, SSE2, SSE3, SSSE3, SSE4.1, SSE4.2; register width: 128 bits; data types: single-precision floating point (32 bits), double-precision floating point (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 type: single-precision floating point (32 bits); Purpose: mainly used for floating-point operations, especially for games and multimedia applications; Registers: MM0 to MM7.

[0013] 3. ARM's SIMD extension: Name: NEON register; Width: 64 bits and 128 bits; Data types: single-precision floating point (32 bits), double-precision floating point (64 bits), integer (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] Take NEON as an example. NEON technology is an advanced SIMD (single instruction multiple data) architecture for ArmCortex-A series processors. It can accelerate multimedia and signal processing algorithms such as video encoders / decoders, 2D / 3D graphics, games, audio and voice processing, image processing, phone calls and sounds. NEON instructions perform "packed SIMD" processing, registers are considered vectors of elements of the same data type, and data types support: signed / unsigned 8-bit, 16-bit, 32-bit, 64-bit.

[0015] There are 32 NEON vector registers in ARMv8, with a length of 128 bits. They can be used to process both scalar operands and vector operands. These registers can store multiple operands of the same data type, such as two 64-bit double type operands or four 32-bit float type operands. During the execution of NEON instructions, multiple operands on the same vector register are processed simultaneously, thereby achieving parallel processing of different data. Using vector registers to accelerate programs requires reading data from memory into them in advance, so 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, we are often faced with solving equations with the same distribution of non-zero elements. For example, in the field of electromagnetic computing, when the model and mesh division remain unchanged, the distribution of non-zero elements in the equations generated at different frequencies is the same. Since there are a large number of complex matrices (dense complex matrices, sparse complex matrices, etc.) in the field of electromagnetic numerical analysis that need to be solved linearly, the current technology requires calling a solution function for each solution, and there is no way to solve two or more linear equations in parallel.

[0017] Based on this, the inventors of the present application designed a parallel LU direct solution 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 to achieve accelerated solution.

[0019] In order to achieve the above purpose, the embodiments of the present application are implemented in the following manner:

[0020] In a first aspect, the present application provides a data processing method for SIMD parallel LU direct solution, comprising: obtaining data to be processed, wherein the data to be processed includes 2 m pattern matrix, m∈Z+ , each pattern matrix has the same non-zero element distribution pattern; if the pattern matrix is ​​a complex matrix, each pattern matrix is ​​processed into 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 the permutation matrix corresponding to the pattern matrix are vectorized to obtain a vectorized one-dimensional array for storage; the vectorized one-dimensional array is read using the SIMD extension component, parallel LU decomposition and direct solution are performed, and the direct solution result is output.

[0021] In combination with the first aspect, in a first possible implementation manner of the first aspect, 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, using the following formula for real number processing:

[0022] Kx=b,

[0023]

[0024] 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 The four sub-matrices corresponding to the real number matrix form 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.

[0025] In combination with the first aspect, in a second possible implementation of the first aspect, the real number matrix corresponding to each pattern matrix is ​​permuted to obtain a corresponding transformation matrix and a corresponding permutation matrix, including: for each real number matrix A corresponding to the pattern matrix j , j∈[1,2 m ]: Create a permutation matrix P j , in A j x=b j Multiply both ends by P j,have:

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

[0027] Among them, the permutation matrix P j Initially, it is the unit matrix, Aj is the real number matrix corresponding to the pattern matrix, x j is the unknown vector, b j is the right-hand vector; for the 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.

[0028] In combination with the second possible implementation manner of the first aspect, in a third possible implementation manner of the first aspect, based on 2 m The transformation matrix and permutation matrix corresponding to each pattern matrix are vectorized to obtain a vectorized one-dimensional array for storage, including: 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 is 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.

[0029] In combination with the third possible implementation of the first aspect, in a fourth possible implementation of the first aspect, a SIMD extension component is used to read the vectorized one-dimensional array, and parallel LU decomposition and direct solution are performed, including: using the SIMD extension component to read the vectorized one-dimensional array A, performing parallel LU decomposition, and obtaining 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 jThe 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 Elements at the same position are combined into a set of vectors corresponding to the position, and a vectorized one-dimensional array U is obtained for storage; the vectorized one-dimensional array L, one-dimensional array U and one-dimensional array B are read using the SIMD extension component to perform parallel direct solution calculations.

[0030] In combination with the fourth possible implementation of the first aspect, in a fifth possible implementation of the first aspect, a vectorized one-dimensional array A is read using a SIMD extension component, and a parallel LU decomposition is performed to obtain a lower triangular matrix L corresponding to each pattern matrix. j , upper triangular matrix U j , including: initialization is a 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, and 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 .

[0031] In combination with the fourth possible implementation of the first aspect, in a sixth possible implementation of the first aspect, a SIMD extension component is used to read the vectorized one-dimensional array L, one-dimensional array U and one-dimensional array V, and a parallel direct solution calculation is performed, including: based on the lower triangular matrix L corresponding to each pattern matrix j , upper triangular matrix U j ,have:

[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, and 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 containing 2 m elements, based on 2 in each set of vectors m elements and each intermediate vector yj Perform the corresponding parallel multiplication and addition calculations to solve each unknown vector x j .

[0034] In a second aspect, the present application provides a data processing device for SIMD parallel LU direct solution, including: a data acquisition unit for acquiring data to be processed, wherein the data to be processed includes 2 m pattern matrix, m∈Z + , each pattern matrix has the same non-zero element distribution pattern; a real numbering unit is used to perform real numbering on each pattern matrix when the pattern matrix is ​​a complex number matrix 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; a vectorization unit is used to perform vectorization on the real number matrix based on 2 m The transformation matrix and the permutation matrix corresponding to the pattern matrix are vectorized to obtain the vectorized one-dimensional array for storage; the parallel solution unit is used to 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.

[0035] In a third aspect, an embodiment of the present application provides a storage medium, which is arranged in an electronic device and includes 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 as described in the first aspect or any one of the possible implementation methods of the first aspect.

[0036] In a fourth aspect, an embodiment of the present application provides 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, and when the program instructions are loaded and executed by the processor, the steps of the data processing method for parallel LU direct solution for SIMD as described in the first aspect or any one of the possible implementation methods 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 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 SIMD extension components to achieve parallel LU direct solution and improve the solution efficiency, it is necessary to consider the transformation of complex matrices and the problem of being unable to perform LU decomposition. m pattern matrix, m∈Z +, each pattern matrix has the same non-zero element distribution pattern. 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. If it is a real number matrix, it does not need to be converted to a real number, thereby solving the problem of converting complex number matrices to a real number. By permuting the real number 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 number matrix is ​​transformed, and the permutation matrix is ​​the matrix that reflects the transformation process of the real number matrix. The permutation operation is based on determining the main element of the matrix. By permuting the real number 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 , 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. This permutation operation is based on columns, which can effectively avoid the situation where the LU decomposition cannot be performed or the calculation is wrong because the matrix is ​​singular or nearly singular. m The transformation matrix and permutation matrix corresponding to each pattern matrix are vectorized and stored in a one-dimensional array. The SIMD extension component is then used to read the one-dimensional array and perform parallel LU decomposition and direct solution to output the direct solution result. This solution can effectively utilize the acceleration performance of the SIMD extension component to achieve parallel linear solution calculation, greatly improve computing efficiency, and shorten solution time.

[0039] 2. 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 is j b j The elements at the same position are combined into a set of vectors corresponding to the position, and the vectorized one-dimensional array B is stored. In this way, the elements at the same position of all pattern matrices can be combined into a vector at that 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 it later, the jump step size can be set to 2 accordingly. mTherefore, when performing calculation operations, SIMD instructions can be used to perform "packed SIMD" processing, so that multiple operands on the same vector register are processed simultaneously, thereby realizing parallel processing of different data. Using vector registers to accelerate programs requires reading data from memory into them in advance. Data alignment can significantly improve the performance of SIMD extension components, effectively increase the calculation speed of parallel processing, and thus improve the efficiency of direct solution.

[0040] 3. In the solution process, through permutation, LU decomposition, and solution calculation, it can be illustrated as this process: 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 the LU decomposition and solution calculation process 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 pattern matrix. j , upper triangular matrix U j ; The lower triangular matrix L 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 obtained for storage. 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 vectorized one-dimensional array L, one-dimensional array U and one-dimensional array B are read using the SIMD extension component to perform parallel direct solution calculations. This method can play a great role in the solution process of many pattern matrices (such as sparse complex matrices and dense complex matrices with the same type). Not only can the SIMD extension component be used to achieve parallel processing, but also the aligned data storage scheme can be used to ensure the continuity of the accessed data in the memory, greatly improving the acceleration performance of the SIMD extension component. For example, a test case is provided in this solution, which is two dense complex matrices of size 942*942. The comparison method is to perform LU direct solution on the two matrices respectively without using the SIMD extension component, and record the solution time. When the LU is directly solved, it takes 750.794s, while this solution (taking NEON parallel direct solution as an example) only takes 141.341s, and the efficiency is improved by 81.17%.

[0041] In order to make the above-mentioned objects, features and advantages of the present application more obvious and easy to understand, preferred embodiments are specifically cited below and described in detail with reference to the attached drawings. BRIEF DESCRIPTION OF THE DRAWINGS

[0042] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the drawings required for use in 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 therefore should not be regarded as limiting the scope. For ordinary technicians in this field, other related drawings can be obtained based on these drawings without paying creative work.

[0043] Figure 1 A flowchart of a data processing method for SIMD-oriented parallel LU direct solution provided in an embodiment of the present application.

[0044] Figure 2 is the pattern matrix K 1 After real numbering, we get the real matrix A 1 Schematic diagram of the process.

[0045] Figure 3 is the pattern matrix K 2 After real numbering, we get the real matrix A 2 Schematic diagram of the process.

[0046] Figure 4 A is a real matrix 1 After the permutation operation, the transformation matrix A is obtained 1 Schematic diagram of the process.

[0047] Figure 5 A is a real matrix 2 After the permutation operation, the transformation matrix A is obtained 2 Schematic diagram of the process.

[0048] Figure 6 To transform the matrix A 1 ′ and transformation matrix A 2 ′ is a schematic diagram for vectorization.

[0049] Figure 7 A schematic diagram of combining vectors into a one-dimensional array for storage. DETAILED DESCRIPTION

[0050] The technical solutions in the embodiments of the present application will be described below in conjunction with the drawings in the embodiments of the present application.

[0051] In the field of numerical analysis (such as electromagnetic field numerical analysis, mechanical numerical analysis, etc.), there is often a need to solve matrices. Among these matrices that need to be solved (there are real number matrices and complex number matrices, such as dense complex number matrices, sparse complex number matrices, etc.), there are many pattern matrices, where the pattern matrix represents a matrix with the same non-zero element distribution pattern. For example, in the field of electromagnetic calculation, when the model and grid division remain unchanged, the distribution position of the non-zero elements of the equations generated at different frequencies is the same. To solve these matrices, direct solution methods, iterative solution methods and other solution methods can be used. However, the existing direct solution method requires calling a solution function every time a solution is obtained, and there is no way to solve two or more linear equations in parallel.

[0052] Based on this, this embodiment proposes a data processing method for SIMD-oriented parallel LU direct solution, using SIMD extension components to achieve the simultaneous horizontal solution of multiple equations with the same non-zero element distribution. For example, a model in the electromagnetic field has been mesh discretized, and the electromagnetic response at two frequencies of 100Hz and 10Hz is solved. Each equation group needs to be solved, and the non-zero element distribution positions of these two equation groups are exactly the same, resulting in two pattern matrices. Therefore, it is suitable to use the method proposed by the present invention to vectorize the solution process horizontally to achieve parallel solution. .

[0053] See also Figure 1 , Figure 1 A flow chart of the data processing method for SIMD-oriented parallel LU direct solution provided in this embodiment.

[0054] First, step S10 may be executed.

[0055] Step S10: Obtain 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.

[0056] In this embodiment, the data to be processed, that is, the pattern matrix to be solved (usually also including boundary conditions, written in the form of Kx=b), can be obtained. The pattern matrix to be processed once needs to satisfy 2 m The number of +, each pattern matrix has the same non-zero element distribution pattern. At present, the latest ARM9 platform supports the execution of 4 double-precision instructions at the same time, and it may be expanded in the future. Although the technical difficulty is high and the double-precision instructions that can be processed simultaneously are limited (currently the level of executing 2-4 double-precision instructions at the same time), it still has great application potential, and the idea of ​​this scheme can theoretically achieve higher-level parallel processing to greatly improve processing efficiency. This embodiment is introduced by taking 2 pattern matrices as an example for convenience of explanation, which should not be regarded as a limitation of this application.

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

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

[0059] Step S20: 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.

[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 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 The four sub-matrices corresponding to the real matrix form are recorded as the first sub-matrix, the second sub-matrix, the third sub-matrix, and the fourth sub-matrix. Then 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 matrix corresponding to the pattern matrix. Assuming K is an n×n complex 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 real matrix A formed j That is, 2n×2n, where A jis the real matrix corresponding to the jth pattern matrix, j∈[1,2 m ].

[0064] For example Figure 2 and Figure 3 As shown, Figure 2 The pattern matrix K on the left 1 After real numbering, the real number matrix A consisting of the first submatrix, the second submatrix, the third submatrix, and the fourth submatrix on the right is obtained. 1 ;Will Figure 3 The pattern matrix K on the left 2 After real numbering, the real number matrix A consisting of the first submatrix, the second submatrix, the third submatrix, and the fourth submatrix on the right is obtained. 2 .

[0065] After obtaining the real number matrix corresponding to each pattern matrix, step S30 may be further performed.

[0066] Step S30: 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.

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

[0068] For example, for each pattern matrix corresponding to the real number matrix A j :

[0069] We can create a permutation matrix P j (Initially the identity matrix), in A j x=b j Multiply both ends by P j ,have:

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

[0071] 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.

[0072] The specific operations are as follows:

[0073] 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 (excluding the elements in the row where the pivot position of the determined pivot is located, that is, assuming that the pivot at (1, 1) has been determined and the pivot at (2, 2) is currently being determined, it can only be determined from the remaining data of the undetermined pivot) in the pivot order, and exchange the row containing the element with the row where the pivot position is located (if there are the same absolute values, exchange them with the row with the smallest row number), 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 ,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′ For the real matrix A j The transformation matrix after the permutation operation, and the permutation matrix P at this time j , is also gradually updated with the transformation operation. No other characters and the initial permutation matrix P are used here. j In addition, in rare cases, the last row or rows may have a pivot of 0 due to limited or no exchange. In this case, the direct solution cannot be used (the system may prompt that the matrix is ​​a singular matrix and LU decomposition cannot be performed), so other methods (such as iterative method) need to be considered for solution.

[0076] like Figure 4 and Figure 5 As shown, the pattern matrix K 1 The corresponding real matrix A 1 After the permutation operation, the corresponding transformation matrix A is obtained 1′ ,like Figure 4 ; The pattern matrix K 2 The corresponding real matrix A 2 After the permutation operation, the corresponding transformation matrix A is obtained 2′ ,like Figure 5 .

[0077] For each pattern matrix corresponding to the real matrix A j After the permutation operation, 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 2 m The transformation matrix and permutation matrix corresponding to the pattern matrix are vectorized to obtain a vectorized one-dimensional array for storage.

[0079] In order to align the data and further improve the acceleration performance of the SIMD extension component, the transformation matrix A corresponding to all pattern matrices can be 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.

[0080] For example, the transformation matrix A corresponding to each pattern matrix can be j′ Elements at the same position are combined into a set of vectors corresponding to that position (including 2 m elements), and then form a one-dimensional array for storage in the order of row and column numbers. This ensures the continuity of the SIMD expansion component when reading and effectively exerts the acceleration performance.

[0081] like Figure 6 As shown, the transformation matrix A 1′ And the transformation matrix A 2′ The elements at the same position in the matrix are combined into a set of vectors corresponding to the position, such as position (1, 1), which is the position of the first row and the first column, and the transformation matrix A 1′ The element at this position is -6, and the transformation matrix A 2′ The element at this position is 3, which is combined into a set of vectors (-6, 3). The same combination is performed on other positions in the first row to determine the corresponding vectors, such as Figure 7 As shown (showing the first row of vectors forming an array), the two elements in the dashed box are a set of vectors (the gray background represents the transformation matrix A 1′ The corresponding position of the value, the white background represents the transformation matrix A 2′ With the transformation matrix A 1′ The values ​​of the aligned positions are used to facilitate the subsequent SIMD extension components to read the information of two transformation matrices at a time for LU decomposition), and the vectorized one-dimensional array A is obtained in order (for example, row 1, column 1 - row 1, column 2... row 2, column 1... row 2n, column 2n) for storage. Use double type to save, 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 SIMD extension components. Here, if 4 transformation matrices are vectorized, a 128-bit vector can contain 4 32-bit float type data to accommodate parallel processing of more matrices.

[0082] Similarly, the right-hand vector P corresponding to all pattern matrices can be 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.

[0083] Thus, the processing of data (2 m After the pre-processing of the pattern matrix), step S50 can be executed.

[0084] Step S50: 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.

[0085] In this embodiment, the SIMD extension component can be used 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 .

[0086] Due to modern processor architectures (mainstream ones include x86 architecture and arm architecture, ARM NEON is the SIMD instruction set under the ARM platform, and the x86 platform has MMX, SSE, AVX and other instruction sets, and some other architectures also have some SIMD instruction sets), it is possible to support the simultaneous execution of multiple instructions through its SIMD extension components (vectorization components). This embodiment takes NEON as an example, which should not be regarded as a limitation of the present application.

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

[0088] Assume that the order of A is the principal subform D i ≠0(i=1:n), then A=LU, note:

[0089]

[0090] By 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 sequence using the above formula.

[0095] This solution is to achieve parallelism, you can initialize is a double precision type, and the value is set to 0.0; each group of vectors in the one-dimensional array A is read one by one through the NEON register, and each group 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 .

[0096] The specific decomposition algorithm is designed (taking the parallel calculation of two matrices as an example), and the code and explanation are as follows:

[0097]

[0098]

[0099]

[0100] This allows parallel processing of LU decomposition to be achieved, and the updated lower triangular matrix L corresponding to each pattern matrix is ​​obtained. j , upper triangular matrix U j .

[0101] Get the updated lower triangular matrix L corresponding to each pattern matrix j , upper triangular matrix U j After that, the lower triangular matrix L corresponding to all pattern matrices can be 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 vectorized storage method is described in detail in the previous text and will not be described here.

[0102] Afterwards, SIMD extension components (such as NEON registers) can be used to read the vectorized one-dimensional array L, one-dimensional array U, and one-dimensional array B, perform parallel direct solution calculations, and output direct solution results.

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

[0104]

[0105] Among them, A j′ x j =Pj b j Refers to the replacement process, (L j U j )x j =P j b j Refers to the LU decomposition process, Refers to the subsequent solution process. In the solution process, you need to first pass L j y j =P j b j Solving for the intermediate vector y j , and then the intermediate vector y j Substitute U j x j =y j , to solve for the unknown vector x j .

[0106] Specifically, each intermediate vector y can be initialized j , j∈[1,2 m ], read each set of vectors in the one-dimensional array L and each set of vectors in the one-dimensional array B one by one through the NEON register, 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 . For permutation, each set of vectors in the one-dimensional array U can be read one by one through the NEON register. 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 .

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

[0108]

[0109]

[0110] Thus, each unknown vector x can be calculated j , that is, the solution corresponding to each pattern matrix. The output solution can also be further converted and restored to complex form. For example, since x is of 2×n dimension and the original complex matrix dimension is 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. According to this rule, the corresponding complex solution can be obtained by performing complex reduction operation.

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

[0112] Table 1. Comparison results of this embodiment and conventional solution

[0113] Unit\Method NEON parallel direct solution LU direct solution NEON efficiency improvement Time(s) 141.341 750.794 81.17％

[0114] It can be seen that LU direct solution requires 750.794s, while this solution (NEON parallel direct solution) only requires 141.341s, and the efficiency is improved by 81.17%.

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

[0116] 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.

[0117] The real number conversion unit is used to perform real number conversion 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.

[0118] The 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.

[0119] 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.

[0120] 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.

[0121] In this embodiment, the real number conversion unit is specifically used to: for each pattern matrix, use the following formula to perform real number conversion:

[0122] Kx=b,

[0123]

[0124] K is the current mode matrix, x is the unknown vector, b is the right end vector, Kr 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 The four sub-matrices corresponding to the real number matrix form 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.

[0125] In this embodiment, the matrix permutation unit is specifically used to: for each real number matrix A corresponding to the pattern matrix j , j∈[1,2 m ]: Create a permutation matrix P j , in A j x=b j Multiply both ends by P j ,have:

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

[0127] Among them, the permutation matrix P j Initially, it is the identity matrix, A j is the real number matrix corresponding to the pattern matrix, j is the unknown vector, b j is the right-hand vector; for the 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.

[0128] In this embodiment, the vectorization unit is specifically used to:m The transformation matrix and permutation matrix corresponding to each pattern matrix are vectorized to obtain a vectorized one-dimensional array for storage, including: 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 is 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.

[0129] In this embodiment, the parallel solving unit is specifically used 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 , 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 Elements at the same position are combined into a set of vectors corresponding to the position, and a vectorized one-dimensional array U is obtained for storage; the vectorized one-dimensional array L, one-dimensional array U and one-dimensional array B are read using the SIMD extension component to perform parallel direct solution calculations.

[0130] In this embodiment, the parallel solution unit is specifically used to: initialize is a 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, and 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 .

[0131] In this embodiment, the parallel solving unit is specifically used to: based on the lower triangular matrix L corresponding to each pattern matrix j , upper triangular matrix U j ,have:

[0132]

[0133] 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, and 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 containing 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 .

[0134] This embodiment further provides a storage medium, which is arranged in an electronic device and includes a stored program. When the program is running, the electronic device where the storage medium is located is controlled to execute the data processing method for SIMD-oriented parallel LU direct solution of this embodiment.

[0135] In addition, this embodiment further provides an electronic device, including 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, and 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 of this embodiment are implemented.

[0136] In summary, the embodiment of the present application provides a data processing method and device for parallel LU direct solution for SIMD. In the field of electromagnetic numerical analysis (or other numerical analysis fields, such as mechanical numerical analysis fields), 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 SIMD extension components to realize parallel LU direct solution and improve the solution efficiency, it is necessary to consider the transformation problem of complex matrices and the problem of being unable to perform LU decomposition. Accordingly, this solution obtains the data to be processed, including 2 m pattern matrix, m∈Z + , each pattern matrix has the same non-zero element distribution pattern. 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. If it is a real number matrix, it does not need to be converted to a real number, thereby solving the problem of converting complex number matrices to a real number. By permuting the real number 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 number matrix is ​​transformed, and the permutation matrix is ​​the matrix that reflects the transformation process of the real number matrix. The permutation operation is based on determining the main element of the matrix. By permuting the real number matrix A jFor 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 , 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. This permutation operation is based on columns, which can effectively avoid the situation where the LU decomposition cannot be performed or the calculation is wrong because the matrix is ​​singular or nearly singular. m The transformation matrix and permutation matrix corresponding to each pattern matrix are vectorized and stored in a one-dimensional array. The SIMD extension component is then used to read the one-dimensional array and perform parallel LU decomposition and direct solution to output the direct solution result. This solution can effectively utilize the acceleration performance of the SIMD extension component to achieve parallel linear solution calculation, greatly improve computing efficiency, and shorten solution time.

[0137] 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 is j b j The elements at the same position are combined into a set of vectors corresponding to the position, and the vectorized one-dimensional array B is stored. In this way, the elements at the same position of all pattern matrices can be combined into a vector at that 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 it later, the jump step size can be set to 2 accordingly. m Therefore, when performing calculation operations, SIMD instructions can be used to perform "packed SIMD" processing, so that multiple operands on the same vector register are processed simultaneously, thereby realizing parallel processing of different data. Using vector registers to accelerate programs requires reading data from memory into them in advance. Data alignment can significantly improve the performance of SIMD extension components, effectively increase the calculation speed of parallel processing, and thus improve the efficiency of direct solution.

[0138] In the solution process, through permutation, LU decomposition, and solution calculation, it can be illustrated as this process: 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 the LU decomposition and solution calculation process 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 pattern matrix. j , upper triangular matrix U j ; The lower triangular matrix L 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 obtained for storage. 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 vectorized one-dimensional array L, one-dimensional array U and one-dimensional array B are read using the SIMD extension component to perform parallel direct solution calculations. This method can play a great role in the solution process of many pattern matrices (such as sparse complex matrices and dense complex matrices with the same type). Not only can the SIMD extension component be used to achieve parallel processing, but also the aligned data storage scheme can be used to ensure the continuity of the accessed data in the memory, greatly improving the acceleration performance. For example, a test case is provided in this solution, which is two dense complex matrices of size 942*942. The comparison method is to perform LU direct solution on the two matrices respectively without using the SIMD extension component, and record the solution time. When the LU is directly solved, it takes 750.794s, while this solution (taking NEON parallel direct solution as an example) only takes 141.341s, and the efficiency is improved by 81.17%.

[0139] The above description is only an embodiment of the present application and is not intended to limit the protection scope of the present application. For those skilled in the art, the present application may have various modifications and variations. Any modification, equivalent replacement, improvement, 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: Initialize L 1 , …, U 1 , …, 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 L 1 , …, U 1 , …, 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