A Parallel Supernode Sparse Triangular Solving Method for Structural Finite Element Stiffness Equations

Abstract

This application belongs to the field of structural finite element analysis, and particularly relates to a parallel supernode sparse triangular solution method for structural finite element stiffness equations. The method includes: performing finite element division on an aircraft structure to obtain the stiffness matrix and load matrix of the finite element model of the aircraft structure, and constructing the structural finite element stiffness equation Ax = B ; performing a first-row permutation operation on the load matrix B ; performing forward substitution based on the supernode data structure; performing backward substitution based on the supernode data structure; performing a second-row permutation operation on the load matrix B that is opposite to the first-row permutation operation to obtain the displacement matrix x ; performing finite element analysis of the aircraft structure according to the displacement matrix x . This application adopts triangular matrix solution and matrix multiply-add operations for a new parallel sparse triangular solution. Taking supernodes as units, forward substitution and backward substitution are performed for each supernode, greatly improving the efficiency of forward substitution and backward substitution.

Description

Technical Field

[0001] This application belongs to the field of structural finite element analysis, and particularly relates to a parallel supernode sparse triangular solution method for structural finite element stiffness equations. Background Art

[0002] The strength design of aircraft structures aims to ensure that the aircraft can withstand the effects of various working environments such as loads, vibrations, and temperatures during flight, avoid structural failures, and maintain safe operation. Structural finite element analysis is widely used in the static strength analysis and dynamic strength analysis of aircraft structures.

[0003] By performing sparse triangular solution on the structural finite element stiffness equation, the structural displacement field can be obtained. This structural displacement field can be used for aircraft structural deformation analysis and aircraft structural dynamics analysis, etc., to support the quantitative evaluation of aircraft structural performance and guide the improvement of aircraft structural design. In the traditional sparse triangular solution process, when using the lower triangular matrix for forward substitution and back substitution, the matrix is usually updated column by column, and the operation efficiency of this method is relatively low.

[0004] Therefore, it is desirable to have a technical solution to overcome or at least mitigate at least one of the above defects of the prior art. Summary of the Invention

[0005] The purpose of this application is to provide a parallel supernode sparse triangular solution method for structural finite element stiffness equations to solve at least one problem existing in the prior art.

[0006] The technical solution of this application is as follows:

[0007] A parallel supernode sparse triangular solution method for structural finite element stiffness equations includes:

[0008] Step 1: Perform finite element division on the aircraft structure, obtain the stiffness matrix and load matrix of the aircraft structure finite element model, and construct the structural finite element stiffness equation Ax = B , where A is the stiffness matrix, x is the displacement matrix, B is the load matrix;

[0009] Step 2: Perform the first row permutation operation on the load matrix B ;

[0010] Step 3: Perform forward substitution based on the supernode data structure;

[0011] Step 4: Perform back substitution based on the supernode data structure;

[0012] Step 5: Perform BPerform a second row permutation operation opposite to the first row permutation operation to obtain a displacement matrix x ;

[0013] Step Six: Perform finite element analysis of the aircraft structure according to the displacement matrix x .

[0014] In at least one embodiment of the present application, in Step One, perform finite element division on the aircraft structure to obtain the stiffness matrix and load matrix of the finite element model of the aircraft structure, including:

[0015] Perform finite element division on the aircraft structure to obtain a finite element model of the aircraft structure;

[0016] Obtain the element stiffness matrix of the finite element model of the aircraft structure, and superimpose all the element stiffness matrices to obtain a stiffness matrix;

[0017] Obtain the external load of the aircraft structure, and map the external load to the finite element nodes of the finite element model of the aircraft structure to obtain a load matrix.

[0018] In at least one embodiment of the present application, the aircraft structure is a fuselage or a wing box section.

[0019] In at least one embodiment of the present application, in Step Two, perform a first row permutation operation on the load matrix B , including:

[0020] Perform a first row permutation operation on the element B in the i th row and j th column of the load matrix b ij . After permutation, it is located in the iperm i th row and j th column, where iperm is the inverse permutation sequence obtained by METIS reordering of the stiffness matrix A ;

[0021] Save the load matrix B in column-major order dense form.

[0022] In at least one embodiment of the present application, in Step Three, perform forward substitution based on the supernode data structure, including:

[0023] Add all the supernodes A in the stiffness matrix J that have 0 descendant supernodes to the task pool and create multiple threads;

[0024] Judge B whether 1 is a zero matrix. If not, for the load matrix​B Perform triangular matrix solution operation:

[0025] ;

[0026] Perform matrix multiply-add operation on the load matrix B :

[0027] ;

[0028] Provide a temporary array for each thread TB , then the transformation is:

[0029] ;

[0030] Perform atomic operation to obtain:

[0031] ;

[0032] Among them, B 1, B 2, B 3, B 4, B 5, B 6 are respectively the segment pointers of the load matrix B , C 1, C 2, C 3, C 4, C 5, C 6 are respectively the segment pointers of the supernode J , and C 1 is the lower triangular matrix.

[0033] In at least one embodiment of the present application, in step three, when performing forward substitution based on the supernode data structure, it further includes:

[0034] Supplement the task pool, and determine whether all descendant supernodes of the parent supernode of the supernode J have completed the forward substitution task. If so, add the parent supernode to the task pool.

[0035] In at least one embodiment of the present application, in step four, when performing backward substitution based on the supernode data structure, it includes:

[0036] Add all supernodes A in the stiffness matrix that have 0 ancestor supernodes J to the task pool, and create multiple threads;

[0037] Perform matrix multiply-add operation on the load matrix B :

[0038] ;

[0039] Perform triangular matrix solution operations on the load matrix B :

[0040] ;

[0041] Among them, B 1. B 2. B 3. B 4. B 5. B 6 are respectively the segment pointers of the load matrix B , C 1. C 2. C 3. C 4. C 5. C 6 are respectively the segment pointers of the supernode J , and C 1 is the lower triangular matrix.

[0042] In at least one embodiment of the present application, in step four, when performing back substitution based on the supernode data structure, it further includes:

[0043] Supplement the task pool, and determine whether all ancestor supernodes of the descendant supernodes of the supernode J have completed the back substitution task. If so, add the descendant supernodes to the task pool.

[0044] In at least one embodiment of the present application, determining whether all ancestor supernodes of the descendant supernodes of the supernode J have completed the back substitution task includes:

[0045] If the number of ancestor supernodes of the descendant supernodes of the supernode J that have completed the back substitution task is equal to the number of all ancestor supernodes of the descendant supernodes, then all ancestor supernodes of the descendant supernodes of the supernode J have completed the back substitution task.

[0046] The invention has at least the following beneficial technical effects:

[0047] The parallel supernode sparse triangular solution method for the structural finite element stiffness equation of the present application adopts triangular matrix solution and matrix multiply-add operations to perform a new parallel sparse triangular solution. Taking supernodes as units, forward substitution and back substitution are performed for each supernode, greatly improving the efficiency of forward substitution and back substitution. BRIEF DESCRIPTION OF THE DRAWINGS

[0048] Figure 1It is a flowchart of the parallel supernode sparse triangular solution method for the structural finite element stiffness equation according to an embodiment of the present application;

[0049] Figure 2 It is a schematic diagram of forward substitution based on the supernode data structure according to an embodiment of the present application;

[0050] Figure 3 It is a schematic diagram of backward substitution based on the supernode data structure according to an embodiment of the present application. Specific embodiments

[0051] To make the purpose, technical solutions, and advantages of the implementation of the present application clearer, the technical solutions in the embodiments of the present application will be described in more detail below with reference to the accompanying drawings in the embodiments of the present application. In the drawings, the same or similar reference numerals denote the same or similar elements or elements with the same or similar functions from beginning to end. The described embodiments are some, but not all, of the embodiments of the present application. The embodiments described below by referring to the drawings are exemplary and are intended to explain the present application and should not be construed as limiting the present application. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present application without creative efforts fall within the scope of protection of the present application. The embodiments of the present application will be described in detail below with reference to the accompanying drawings.

[0052] In the description of the present application, it should be understood that the orientation or positional relationship indicated by the terms "center", "longitudinal", "transverse", "front", "rear", "left", "right", "vertical", "horizontal", "top", "bottom", "inner", "outer", etc. is based on the orientation or positional relationship shown in the drawings, and is only for the convenience of describing the present application and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and thus should not be construed as limiting the scope of protection of the present application.

[0053] The following combines the attached Figures 1 to 3 to further elaborate on the present application.

[0054] The present application provides a parallel supernode sparse triangular solution method for the structural finite element stiffness equation, as Figure 1 shown, including the following steps:

[0055] Step 1: Perform finite element division on the aircraft structure, obtain the stiffness matrix and load matrix of the aircraft structure finite element model, and construct the structural finite element stiffness equation Ax = B , where A is the stiffness matrix, x is the displacement matrix, B is the load matrix;

[0056] Step 2: Perform the first row permutation operation on the load matrix B ;

[0057] Step 3: Perform forward substitution based on the supernode data structure

[0058] Step 4: Perform backward substitution based on the supernode data structure

[0059] Step 5: Perform the second row permutation operation on the load matrix B which is the opposite of the first row permutation operation to obtain the displacement matrix x ;

[0060] Step 6: Perform finite element analysis of the aircraft structure according to the displacement matrix x ;

[0061] In the parallel supernode sparse triangular solution method for the structural finite element stiffness equation of the present application, in Step 1, the aircraft structure is divided into finite elements to obtain the stiffness matrix and load matrix of the aircraft structure finite element model, which specifically includes:

[0062] Divide the aircraft structure into finite elements to obtain the aircraft structure finite element model

[0063] Obtain the element stiffness matrix of the aircraft structure finite element model, and superimpose all the element stiffness matrices to obtain the stiffness matrix

[0064] Obtain the external load of the aircraft structure, and map the external load to the finite element nodes of the aircraft structure finite element model to obtain the load matrix

[0065] In this embodiment, the aircraft structure is the fuselage or wing box section

[0066] In the parallel supernode sparse triangular solution method for the structural finite element stiffness equation of the present application, for the structural finite element stiffness equation Ax = B , according to the lower triangular matrix and upper triangular matrix obtained from the numerical decomposition in the Cholesky decomposition process, perform two substitution solution processes respectively from front to back forward and from back to front backward, and the equation solution can be obtained. For the Cholesky decomposition, due to symmetry, only the lower triangular matrix needs to be decomposed, and both the forward substitution and backward substitution are performed on this lower triangular matrix. Before performing the substitution solution, the numerical decomposition has been completed, and the data is saved in the form of block supernodes

[0067] In the parallel supernode sparse triangular solution method for the structural finite element stiffness equation of the present application, since the stiffness matrix A has undergone a symmetric permutation operation after METIS reordering, in order to obtain the correct equation solution, the right-hand side load matrix B also needs to be performed the same as the stiffness matrixA The same row permutation operation.

[0068] Specifically, in step two, for the load matrix B perform the first row permutation operation, including:

[0069] For the load matrix B in the i row and the j column elements b ij perform the first row permutation operation. After permutation, it is located in the iperm i row and the j column, where iperm is the inverse permutation sequence obtained by METIS reordering of the stiffness matrix A ;

[0070] Save the load matrix B in column-major order dense form.

[0071] For the parallel supernode sparse triangular solution method of the structural finite element stiffness equation in this application, to make full use of the characteristics of supernodes, a parallel forward substitution algorithm and a parallel backward substitution algorithm based on dense matrix operations are proposed. The core operations are completed by the high-performance multi-threaded level-3 BLAS function and the subroutines dgemm() and dtrsm() of the LAPACK function, and the task pool uses the OpenMP task pool.

[0072] In the preferred embodiment of this application, in step three, perform forward substitution based on the supernode data structure, including:

[0073] Add all supernodes in the stiffness matrix A that have 0 descendant supernodes J to the task pool and create multiple threads;

[0074] Judge B whether 1 is a zero matrix. If not, perform triangular matrix solution operation on the load matrix B :

[0075] (1)

[0076] Perform matrix multiply-add operation on the load matrix B :

[0077] (2)

[0078] Provide a temporary array TB for each thread, then it is transformed into:

[0079] (3)​

[0080] The atomic operation results in:

[0081] (4)

[0082] Among them, B 1, B 2, B 3, B 4, B 5, B 6 are respectively the segment pointers of the load matrix B , C 1, C 2, C 3, C 4, C 5, C 6 are respectively the segment pointers of the supernode J , and C 1 is a lower triangular matrix.

[0083] As Figure 2 shown, in the initial stage, any supernode without descendant supernodes can start forward substitution. When performing the forward substitution corresponding to a certain supernode, the segment pointer C 1 corresponding to the segment pointer B 1 part has no possibility of multiple threads operating simultaneously, so there will be no data competition in the triangular matrix solution operation. However, since the load matrix B parts corresponding to the non-zero elements below the diagonal of different supernodes may overlap, data competition may occur. At this time, a temporary array TB needs to be provided for each thread, and its original data is initialized to 0, and then equation (2) is changed to equation (3). After the calculation is completed, atomic operations can be used to avoid data competition.

[0084] In this embodiment, the parallel forward substitution algorithm is specifically as follows:

[0085]

[0086] When the forward substitution task for the supernode J is completed, the task pool is replenished. When adding tasks to the task pool, it is necessary to judge whether all the descendant supernodes of the parent supernode J of the supernode P J have completed the corresponding forward substitution tasks. If so, the parent supernode P J is added to the task pool. In addition, since the right-hand side load matrix B usually also has sparsity, before performing the triangular matrix solution and matrix multiply-add operations, it is necessary to judge the segment pointer​​B Whether 1 is a zero matrix. If so, there is no need to perform the subsequent operations.

[0087] In a preferred embodiment of the present application, the backward substitution is similar to the forward substitution, except that the update order is from N column - 1 to column 0, and it is performed in reverse from back to front. Specifically, in step four, the backward substitution based on the super - node data structure is performed, including:

[0088] Add all the super - nodes in the stiffness matrix A that have 0 ancestor super - nodes J to the task pool and create multiple threads;

[0089] Perform matrix multiply - add operations on the load matrix B :

[0090] (5)

[0091] Perform triangular matrix solution operations on the load matrix B :

[0092] (6)

[0093] Among them, B 1, B 2, B 3, B 4, B 5, B 6 are respectively the segment pointers of the load matrix B , C 1, C 2, C 3, C 4, C 5, C 6 are respectively the segment pointers of the super - node J , and C 1 is a lower triangular matrix.

[0094] For the backward substitution of super - node J, first perform matrix multiply - add operations, and then perform a triangular matrix solution operation once.

[0095] As Figure 3 shown, in the initial stage, all the super - nodes in the stiffness matrix A that have no ancestor super - nodes J should be added to the task pool. Since only the segment pointer corresponding to the super - node segment pointer C 1 is changed in each operation B 1, there will be no data competition problems in both matrix multiply - add and triangular matrix solution operations.

[0096] In this embodiment, the parallel reverse substitution algorithm is specifically as follows:

[0097]

[0098] After the reverse substitution task for the supernode J is completed, the task pool is replenished. When adding tasks to the task pool, if the supernode J has SJ_N J descendant supernodes, for the descendant supernodes J of the supernode SJ J k , only when all the ancestor supernodes of the descendant supernode SJ J k have completed the reverse substitution task, can it be put into the task pool. Among them, judging whether all the ancestor supernodes of the descendant supernode J of the supernode SJ J k have completed the reverse substitution task includes: if the number of ancestor supernodes of the descendant supernode J of the supernode SJ J k that have completed the reverse substitution task is equal to the number of all ancestor supernodes of the descendant supernode SJ J k , then all the ancestor supernodes of the descendant supernode J of the supernode SJ J k have completed the reverse substitution task.

[0099] In the parallel supernode sparse triangular solution method for the structural finite element stiffness equation of this application, after the forward substitution and reverse substitution are completed, at this time, the equation solution is already included in the memory space of the load matrix B . If the final equation solution is to be obtained, only a row permutation operation opposite to the row permutation operation in step two needs to be performed on the load matrix B .

[0100] In the parallel supernode sparse triangular solution method for the structural finite element stiffness equation of this application, finally, in step six, according to the obtained equation solution of the structural finite element stiffness equation, that is, the displacement matrix x , aircraft structural finite element analysis is carried out, so as to realize the evaluation of aircraft structural performance, support multiple aspects such as aircraft structural design optimization, safety verification, and maintenance plan formulation, which is crucial for ensuring the safe, efficient, and economic operation of the aircraft. ​​​​​​​​​​​​​

[0101] The parallel supernode sparse triangular solution method for the structural finite element stiffness equation of the present application is based on the supernode data structure of the stiffness matrix A When performing forward substitution and backward substitution, each time the supernode obtained by symbolic factorization is used as a unit, the load matrix is updated supernode by supernode B A new parallel sparse triangular solution is carried out by using triangular matrix solution and matrix multiply-add operations. The dense matrix operations are performed by using the subroutines dgemm() and dtrsm() of the high-performance multithreaded level-3 BLAS function and LAPACK function, which fully exploits the computing potential of the computer cache. In addition, combined with the OpenMP task pool, the parallel potential of the calculation is maximally exerted, and it has excellent calculation efficiency.

[0102] As mentioned above, it is only the specific implementation manner of the present application, but the protection scope of the present application is not limited thereto. Any changes or substitutions that can be easily thought of by those skilled in the art within the technical scope disclosed in the present application should be covered within the protection scope of the present application. Therefore, the protection scope of the present application shall be subject to the protection scope of the claimed rights.

Claims

1. A parallel super-node sparse triangular solution method for structural finite element stiffness equations, characterized in that: include: Step 1: Perform finite element division on the aircraft structure, obtain the stiffness matrix and load matrix of the aircraft structure finite element model, and construct the structural finite element stiffness equation Ax = B ,in, A is the stiffness matrix, x is the displacement matrix, B is the load matrix; Step 2: Load Matrix B Perform the first row permutation operation; Step 3: Perform forward substitution based on the supernode data structure, including: The stiffness matrix A All supernodes in the J Add to the task pool and create multiple threads; judge B 1 is a 0 matrix, if not, for the load matrix B Perform a triangular matrix solve operation: Load matrix B Perform matrix multiplication and addition operations: Provide a temporary array for each thread TB , then it is transformed into: Perform atomic operations to obtain: in, B 1. B 2. B 3. B 4. B 5. B 6 are the load matrices in a thread respectively B The fragment pointer, C 1. C 2. C 3. C 4. C 5. C 6 are super nodes J The fragment pointer of C 1 is a lower triangular matrix; Step 4: Perform reverse substitution based on the supernode data structure, including: The stiffness matrix A All supernodes in that have 0 ancestor supernodes J Add to the task pool and create multiple threads; Load matrix B Perform matrix multiplication and addition operations: Load matrix B Perform a triangular matrix solve operation: in, B 1. B 2. B 3. B 4. B 5. B 6 are the load matrices in a thread respectively B The fragment pointer, C 1. C 2. C 3. C 4. C 5. C 6 are super nodes J The fragment pointer of C 1 is a lower triangular matrix; Step 5: Load Matrix B Perform a second row permutation operation opposite to the first row permutation operation to obtain the displacement matrix x ; Step 6: According to the displacement matrix x Conduct finite element analysis of aircraft structures.

2. The parallel super-node sparse triangle solution method for the structural finite element stiffness equation according to claim 1 is characterized in that: In step 1, the aircraft structure is divided into finite elements to obtain the stiffness matrix and load matrix of the aircraft structure finite element model, including: Perform finite element division on the aircraft structure to obtain a finite element model of the aircraft structure; Obtaining a unit stiffness matrix of the aircraft structure finite element model, and superimposing all the unit stiffness matrices to obtain a stiffness matrix; The external load of the aircraft structure is acquired, and the external load is mapped to the finite element nodes of the finite element model of the aircraft structure to obtain a load matrix.

3. The parallel super-node sparse triangle solution method for the structural finite element stiffness equation according to claim 2 is characterized in that: The aircraft structure is a fuselage or a wing box section.

4. The parallel super-node sparse triangle solution method for the structural finite element stiffness equation according to claim 1 is characterized in that: In step 2, the load matrix B Perform the first row permutation operation, including: Load matrix B The i Row, No. j Column Elements b ij Perform the first row replacement operation, and the row after replacement is located at iperm [ i ] line, j Column, where iperm is the stiffness matrix A The reverse permuted sequence obtained by METIS reordering; The loading matrix B Stored in dense column-major order.

5. The parallel super-node sparse triangle solution method for the structural finite element stiffness equation according to claim 4 is characterized in that: In step 3, a forward substitution based on the supernode data structure is performed, which also includes: Supplement the task pool and determine the super node J Whether all descendant supernodes of the parent supernode have completed the forward substitution task, if so, the parent supernode is added to the task pool.

6. The parallel super-node sparse triangle solution method for the structural finite element stiffness equation according to claim 5 is characterized in that: In step 4, reverse substitution based on the supernode data structure is performed, which also includes: Supplement the task pool and determine the super node J Whether all ancestor super nodes of the descendant super node have completed the reverse replacement task, if so, the descendant super node is added to the task pool.

7. The parallel super-node sparse triangle solution method for the structural finite element stiffness equation according to claim 6 is characterized in that: Determine super node J Whether all ancestor supernodes of the descendant supernode have completed the reverse replacement task, including: If the super node J The number of ancestor supernodes whose descendant supernode completes the reverse substitution task is equal to the number of all ancestor supernodes of the descendant supernode. J All ancestor supernodes of the descendant supernodes complete the reverse substitution task.

Images