Heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber

Abstract

The invention relates to the technical field of structural material high-performance calculation and numerical simulation, and discloses a metal rubber random vibration damage and fatigue life prediction heterogeneous parallel method which comprises the following steps: inputting finite element model data of a large-scale metal rubber component; a stiffness matrix, a mass matrix and the like are subjected to accelerated calculation by adopting a CPU multi-thread parallel technology, a frequency response function is calculated by combining a multi-process and multi-thread parallel technology, and random vibration response and component fatigue damage based on a frequency domain method are efficiently calculated by utilizing a GPU parallel computing technology. According to the method, a random vibration theory, a frequency domain fatigue analysis method and a finite element calculation method are combined, and meanwhile, multi-process, multi-thread and GPU heterogeneous parallel calculation frameworks are integrated, so that the calculation efficiency is remarkably improved, efficient analysis of random vibration response and fatigue life prediction of the vibration reduction component made of the complex structural material is realized, and the method is suitable for large-scale popularization and application. And an efficient strategy is provided for design optimization and safe service in the fields of aerospace, high-end equipment and the like.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of high-performance computing and numerical simulation of structural materials, in particular to a heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber. BACKGROUND

[0002] Metal rubber is a kind of porous high-energy dissipation structural material, which has excellent properties such as high and low temperature resistance, large load capacity, radiation resistance and performance designability. It is often used as the preferred material for vibration reduction components in aerospace, high-performance weapon equipment and other fields, and metal rubber exhibits significant mechanical performance advantages in extreme environments. However, its service environment usually involves complex random vibration excitation, and the internal metal wires show high nonlinearity due to large-scale friction and wear and fatigue fracture, which brings great challenges to the numerical simulation of vibration response and fatigue damage, and hinders the precise design, fatigue life prediction and engineering safety service evaluation of metal rubber components.

[0003] In the random vibration analysis of metal rubber components, the finite element model often contains a large number of degrees of freedom and complex nonlinear behavior, generating a large amount of calculation data. The traditional serial computing method is low in computing efficiency and cannot meet the real-time requirements, especially when dealing with large-scale unit node stress, power spectral density and fatigue damage calculation, which takes too long time, limiting the engineering application of high-precision simulation. In addition, random vibration fatigue analysis needs to combine with frequency domain method theory, involving spectral moment calculation and multiple damage models, further increasing the computational complexity. Therefore, developing an efficient parallel computing method has become an urgent need to solve the large-scale random vibration fatigue problem of metal rubber and promote its on-demand design and performance optimization.

[0004] Therefore, a heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber is established, which integrates finite element numerical calculation, frequency domain fatigue theory and CPU-GPU collaborative computing framework, and takes large-scale calculation as the starting point, parallelizes eigenvalue solving, stress response power spectral density calculation and fatigue damage assessment, which can significantly improve the computing efficiency, efficiently solve the large-scale random vibration fatigue problem of metal rubber components, and greatly promote the mechanical performance optimization of metal rubber materials, and provide a tool for metal rubber cross-scale performance evolution prediction, on-demand design and engineering safety service evaluation application. SUMMARY

[0005] To overcome the shortcomings of the prior art, the present application provides a heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber, which solves the problem of long time consumption in processing large-scale unit node stress, power spectral density and fatigue damage calculation, which limits the engineering application of high-precision simulation. In addition, random vibration fatigue analysis needs to combine with frequency domain method theory, involving spectral moment calculation and multiple damage models, further increasing the computational complexity.

[0006] To achieve the above object, the application is implemented by the following technical solutions: a metal rubber random vibration damage and fatigue life prediction heterogeneous parallel method, characterized by comprising the following steps:

[0007] S1, inputting finite element model data of a large-scale metal rubber component;

[0008] S2, calculating total stiffness matrix and total mass matrix of the finite element model by using a multi-thread parallel technology based on CPU;

[0009] S3, calculating frequency response function of the finite element model by combining MPI multi-process parallel technology and OpenMP multi-thread parallel technology;

[0010] S4, calculating random vibration response based on the frequency response function by using GPU parallel computing technology, including element node Mises stress power spectral density and node average Mises stress power spectral density;

[0011] S5, calculating random vibration fatigue damage and logarithmic life of element nodes and global nodes by using a variety of frequency domain damage models based on frequency domain fatigue analysis method through GPU parallel computing technology.

[0012] Preferably, the S2 specifically comprises:

[0013] S21, allocating independent memory space for each CPU thread through a multi-thread parallel framework, and allocating elements in the finite element grid to multiple threads for parallel calculation through a dynamic scheduling manner;

[0014] S22, based on the allocated elements, each thread calculates element stiffness matrix and element mass matrix of the element through OpenMP multi-thread parallel calculation in the same parallel region;

[0015] S23, combining local Triplet lists of each thread through a critical section into a global Triplet list, through which element matrices are assembled into corresponding positions of total stiffness matrix and total mass matrix;

[0016] S24, applying boundary conditions by using a large number method, and updating total stiffness matrix and total mass matrix.

[0017] Preferably, in the S22, the calculation manner of the element stiffness matrix after B-Bar correction is:

[0018] ;

[0019] The element mass matrix adopts a lumped mass matrix form, and is calculated through the following formula:

[0020] ;

[0021] wherein, is the element stiffness matrix, is the element mass matrix, is the B-Bar modified geometry matrix, is the material matrix, is the element shape function matrix.

[0022] Preferably, the S3 specifically comprises:

[0023] S31, initializing the MPI distributed computing environment, distributing the total stiffness matrix and the total mass matrix to the remaining processes by the MPI multi-process broadcast of the master process;

[0024] S32, calling the MUMPS open source distributed solver to perform LU decomposition on the total stiffness matrix;

[0025] S33, calling the PARPACK open source library to cooperatively calculate the eigenvalues and eigenvectors of the system by the shift-invert method;

[0026] S34, calculating the displacement response power spectral density of the nodes.

[0027] Preferably, the shift-invert method is to transform the generalized eigenvalue problem into , wherein , to transform the smallest eigenvalue into the largest eigenvalue for fast solving.

[0028] Preferably, the S4 specifically comprises:

[0029] S41, preprocessing the input data and copying from the host memory to the GPU memory;

[0030] S42, configuring the CUDA thread block and grid, each thread processing one element, local node and frequency point;

[0031] S43, calculating the strain using the B-Bar modified B matrix, then calculating the stress at the Gauss integration point through the material D matrix, and combining the shape function matrix to obtain the element node stress;

[0032] S44, calculating the Mises stress response power spectral density of the element node using the Segalman formula;

[0033] S45, calculating the average Mises stress power spectral density of the nodes in parallel through CUDA;

[0034] S46, copying the GPU calculation result from the GPU memory back to the host memory.

[0035] Preferably, the S5 specifically comprises:

[0036] S51, calculating the spectral moment by numerical integration, and further calculating the zero-crossing rate, peak-crossing rate and irregular factor;

[0037] S52, calculating the fatigue damage of the unit node by a plurality of frequency domain fatigue models;

[0038] S53, calculating the average damage and life of the node by GPU parallel computing, wherein the average damage of the node is calculated by using the arithmetic mean;

[0039] S54, outputting the fatigue damage and life calculation result.

[0040] Preferably, the plurality of frequency domain fatigue models comprises the Bendat method, the Wirsching-Light method, the Dirlik method, the Tovo-Benasciutti method, the Alpha=0.75 method and the Zhao-Bake method.

[0041] Preferably, in the S4 and the S5, the GPU parallel computing technology is used and the GPU parallel computing technology is used, both of which comprise the steps of copying the node displacement frequency response function, B-Bar matrix and other data to the GPU memory after flattening the data into continuous arrays from the host memory, and configuring the CUDA thread block and grid.

[0042] Preferably, the formula for calculating the logarithmic life is:

[0043] ;

[0044] wherein, is the logarithmic life, is the damage factor.

[0045] The present application provides a metal rubber random vibration damage and fatigue life prediction heterogeneous parallel method. It has the following beneficial effects:

[0046] The application provides a metal rubber random vibration fatigue damage analysis method based on heterogeneous parallel computing, significantly improves the efficiency and precision of large-scale finite element simulation, and provides an efficient tool for vibration reduction design and safe service of metal rubber components in the fields of aviation and aerospace. First, the calculation of a total stiffness matrix and a total mass matrix is accelerated through multi-threading through the OpenMP technology; then, the parallelization of frequency response calculation is realized by using regional decomposition combined with the MPI and OpenMP technology, and the efficiency bottleneck of traditional serial calculation in a large-scale degree of freedom model is broken. Finally, the calculation of node stress, Mises stress response power spectral density and various frequency domain fatigue damage models such as the Bendat method, the Wirsching-Light method and the Dirlik method is efficiently realized based on the GPU parallel framework, and the calculation of random vibration response and fatigue life prediction is significantly accelerated. The heterogeneous parallel strategy solves the problem of limited computing resources in large-scale random vibration analysis, and provides an efficient and reliable numerical tool for fatigue life prediction, on-demand design and engineering safety service evaluation of metal rubber components. BRIEF DESCRIPTION OF DRAWINGS

[0047] Figure 1 is a flow chart of a metal rubber random vibration damage and fatigue life prediction heterogeneous parallel method of the application;

[0048] Figure 2 is a schematic diagram of an eight-node hexahedral element used in the application;

[0049] Figure 3 is a schematic diagram of a bowl-shaped metal rubber finite element model and load action used in the application;

[0050] Figure 4 is an MPI collection and distribution schematic diagram of the application;

[0051] Figure 5 is a frequency response parallel calculation flow chart of the application;

[0052] Figure 6 is a displacement frequency response function curve obtained by calculating a bowl-shaped metal rubber finite element model used in the application;

[0053] Figure 7 is a Mises stress response power spectral density curve obtained by calculating a bowl-shaped metal rubber finite element model used in the application;

[0054] Figure 8 is a Mises stress response power spectral density cloud chart obtained by calculating a bowl-shaped metal rubber finite element model used in the application;

[0055] Figure 9Damage result and log life result cloud diagram calculated by the bowl-shaped metal rubber finite element model adopted by the present application;

[0056] Figure 10 Is the heterogeneous parallel acceleration ratio result diagram of the present application;

[0057] Figure 11 Is the time-consuming comparison diagram of the present application under the same model scale using serial and parallel computing. DETAILED DESCRIPTION

[0058] The technical solutions in the embodiments of the present application will be clearly and completely described below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor fall within the scope of protection of the present application.

[0059] Embodiment:

[0060] Please refer to the accompanying Figure 1 , the embodiment of the present application provides a metal rubber random vibration damage and fatigue life prediction heterogeneous parallel method, comprising:

[0061] S1, input the finite element model data of any large-scale metal rubber component, the model data includes the geometric information and material properties of eight-node hexahedral elements, for random vibration fatigue analysis.

[0062] Specifically includes the following steps:

[0063] S11, input eight-node hexahedral grid data, including node coordinates and element node connection information, for defining the geometric structure of the finite element model, eight-node hexahedral elements;

[0064] S12, input material property data, including Young's modulus, Poisson's ratio, fatigue strength coefficient, fatigue index and density, for describing the mechanical and fatigue properties of metal rubber;

[0065] S13, input boundary information and other simulation analysis required parameters, including fixed boundary node number, damping ratio, modal order, acceleration power spectral density and excitation direction, for setting random vibration and fatigue analysis, metal rubber finite element simulation model and its boundary condition application method see Figure 3 .

[0066] S2, using multi-thread parallel technology based on CPU, calculate the total stiffness matrix and total mass matrix of the finite element model;

[0067] Specifically includes the following steps:

[0068] S21, each CPU thread is allocated an independent memory space through a multi-thread parallel framework (OpenMP parallel domain, #pragma omp parallel) for storing data such as shape functions, B matrices and B-Bar correction matrices, and then the elements in the finite element grid are allocated to multiple threads for parallel calculation through dynamic scheduling;

[0069] For an eight-node hexahedral element, based on the physical and geometric relationship, the element geometric matrix and shape function matrix are defined, which need to satisfy the following formula:

[0070] ;

[0071] wherein, represents each physical domain after discretization, is the displacement of each node on a single element. is the element geometric matrix, is the element shape function matrix.

[0072] The specific shape function calculation formula is:

[0073] ;

[0074] wherein, is the shape function, , , is the local coordinate value of the eight nodes.

[0075] Then, B-Bar correction is performed, and the geometric matrix is decomposed into a bulk strain part and a deviatoric strain part, Then, the B-Bar correction is performed on the bulk strain part , and the corrected is denoted as , and the specific correction method is as follows:

[0076] Since the strain relationship geometric matrix of the C3D8 element can be written as:

[0077] ;

[0078] For the matrix of a certain node, it can be written as:

[0079] ;

[0080] wherein, , represents the shape function associated with the node , denotes the degree of freedom number, the body strain part may be expressed as:

[0081] ;

[0082] deviatoric strain :

[0083] ;

[0084] For incompressible materials, the body strain part needs to be changed to :

[0085] ;

[0086] The value of can be the value at the centroid of the element, that is:

[0087] ;

[0088] Therefore, The matrix finally becomes , where The average value within the element volume is adopted .

[0089] S22, Since the stiffness matrix of each element is only related to its own node coordinates and material properties, each element is independent of each other, and the calculation has independence, so the element stiffness matrix and the element mass matrix can be calculated in parallel through OpenMP.

[0090] In the same parallel region, each thread calculates the element stiffness matrix and the element mass matrix of the element based on the assigned element through OpenMP multithreading parallel calculation, and the calculation method for the element stiffness matrix after B-Bar correction is:

[0091] ;

[0092] The element mass matrix adopts the form of lumped mass matrix, and is calculated through the following formula:

[0093] ;

[0094] The element stiffness matrix and the mass matrix calculated by each thread are stored in the form of sparse matrix, and the non-zero elements with absolute value greater than a preset threshold are extracted through the Triplet structure to reduce the calculation amount of subsequent assembly. Each thread maps the non-zero elements of the local matrix to the global degrees of freedom according to the degree of freedom index of the element, and generates a local Triplet list;

[0095] S23, merge the local triplet lists of each thread through the critical section into the global triplet list, through which the element matrix is assembled into the corresponding positions of the total stiffness matrix and the total mass matrix, using sparse matrix storage form;

[0096] S24, apply boundary conditions using the large number method, generate stiffness large number coefficients and mass large number coefficients according to the average values of the diagonal elements of the total stiffness matrix and the total mass matrix, and update the diagonal elements of the stiffness matrix and the mass matrix in turn according to the input boundary freedom index and the corresponding boundary value, complete the application of displacement boundary conditions.

[0097] S3, combine MPI multi-process parallel technology and OpenMP multi-thread parallel technology to calculate the frequency response function of the finite element model, including characteristic frequency, characteristic vector and node displacement response;

[0098] Specifically includes the following steps:

[0099] S31, initialize the MPI distributed computing environment, get the number of processes and the current process number, and dynamically adjust the number of OpenMP threads for each MPI process to optimize the allocation of computing resources.

[0100] S32, compress the total stiffness matrix and the total mass matrix by the master process, calculate the number of non-zero elements, and distribute the non-zero values, row indexes and column start indexes of the matrix to the remaining processes through MPI multi-process broadcast. Each process pre-allocates sparse matrix storage space according to the number of non-zero elements broadcasted, synchronizes all process data through MPI barrier, and uses MPI collection and distribution.

[0101] S33, call the MUMPS open source distributed solver to perform LU decomposition on the total stiffness matrix:

[0102] For the generalized eigenvalue problem , let , σ is always 0, and is decomposed into LU, that is, the total stiffness matrix is decomposed into LU.

[0103] Convert the total stiffness matrix into CSR format containing row index, column index and non-zero value through the master process, and broadcast it to other processes. Initialize the MUMPS open source distributed solver and set parameters. All processes call the MUMPS open source distributed solver to decompose the total stiffness matrix into a form suitable for subsequent solution.

[0104] S34, call the PARPACK open source library to solve the characteristic equation of the system by the shift inversion method to obtain the characteristic value and the characteristic vector;

[0105] The free vibration equation of a multi-degree-of-freedom system can be expressed as follows:

[0106] ;

[0107] where, , and represent the mass matrix, the damping matrix and the stiffness matrix of the system respectively, represents the displacement column vector, represents the velocity column vector, represents the acceleration column vector;

[0108] For a multi-degree-of-freedom system without considering damping, the free vibration equation of the system can be expressed as follows:

[0109] ;

[0110] Generally, the free vibration of the system can be regarded as a simple harmonic vibration, so the solution of the above equation is assumed as:

[0111] ;

[0112] where, is the displacement column vector of the nth mode, is the displacement column vector of the nth mode, is the amplitude column vector of the nth mode. , are the natural circular frequency and the initial phase angle of the nth mode respectively; Thus, the mode equation of the system can be obtained as follows:

[0113]

[0114] ;

[0115] For a vibration system, the sufficient and necessary condition for the non-zero solution of the amplitude is:

[0116] ;

[0117] The above equation is the characteristic equation of the system, and the characteristic values and the characteristic vectors of the system can be obtained by solving the equation;

[0118] The parameters of PARPACK (including the number of target characteristic values, the convergence tolerance, the maximum number of iterations and the dimension of Krylov subspace) are initialized, and since PARPACK is better at finding the largest characteristic value, the shift-inversion method is used to call PARPACK;

[0119] The pdsaupd_function of PARPACK is called by each process to cooperatively calculate the characteristic values and the characteristic vectors by the shift-inversion method:

[0120] For generalized characteristic values​​ Problem, can prove Where And in order to find the eigenvalue closest to zero, when Transformed eigenvalues satisfy So the smallest eigenvalue Become the largest eigenvalue At this time, the minimum eigenvalue of the required order can be quickly obtained;

[0121] Each process is calculated by OpenMP parallel computing quality matrix and vector product, as follows:

[0122] ;

[0123] ;

[0124] ;

[0125] The results are broadcast by MPI to maintain consistency between processes, and the pdseupd_function of PARPACK is called to extract the eigenvalues and eigenvectors of the required order. The master process calculates the quality matrix norm of the eigenvector and broadcasts the results of the normalization processing.

[0126] For generalized eigenvalue Problem, Is a symmetric matrix (stiffness matrix), Is a symmetric positive definite matrix (mass matrix), Is the eigenvector, Is the eigenvalue. The norm is defined as The normalized eigenvector x satisfies: .

[0127] S35, according to the specified number of interpolation points, the following interpolation formula is used to generate a frequency point sequence in the characteristic frequency range, as follows:

[0128] ;

[0129] Where, Is the interpolation of the first The frequency, Is the lower limit of the frequency interval, Is the upper limit of the frequency interval, Is the total number of frequency points.

[0130] S36, according to the excitation direction, calculate the excitation vector , and then calculate the first Order modal participation factor by the following formula:

[0131] ;

[0132] The motion differential equation of the multi-degree-of-freedom system under the action of random dynamic load can be expressed as:

[0133] ;

[0134] wherein, , and respectively represent the mass matrix, the damping matrix and the stiffness matrix of the system, therefore, the frequency response function is calculated by the following formula:

[0135] ;

[0136] wherein, is the mass matrix, is the damping matrix, is the stiffness matrix.

[0137] S37, calculate the displacement, velocity and acceleration response power spectral density of the node:

[0138] ;

[0139] ;

[0140] ;

[0141] wherein, is the displacement response power spectral density, is the frequency response function, is the input acceleration power spectral density, is the velocity response power spectral density, is the acceleration response power spectral density;

[0142] The maximum relative error of the displacement frequency response curve calculated by the bowl-shaped metal rubber finite element model and the commercial software is not more than .

[0143] S4, calculate the random vibration response based on the frequency response function by using the GPU parallel computing technology, including the element node stress, the element node Mises stress power spectral density and the node average Mises stress power spectral density;

[0144] Specifically comprising the following steps

[0145] S41, according to the input power spectral density and the interpolation frequency, the weight of each interpolation frequency point is calculated by using the logarithmic interpolation method, so as to provide parameters for subsequent power spectral density calculation;

[0146] S42, preprocessing the input data (including node displacement frequency response function, B-Bar matrix, material D matrix, node coordinates and element node index), flattening to continuous array for GPU memory transmission, and initializing output array, allocating GPU memory for input data and output data (including Gauss integral point stress, element node stress, Mises stress and node average Mises stress), copying flattened data from host memory to GPU memory;

[0147] S43, configuring CUDA thread block and grid, each thread processing one element, local node and frequency point, limiting thread block size to no more than 1024, ensuring composite GPU hardware requirements, and each thread in the thread block processing corresponding node data;

[0148] S44, using B-Bar corrected matrix to calculate strain , and then through material matrix to calculate Gauss integral point stress , combined with shape function matrix to calculate element node stress, and finally write Gauss integral point stress and element node stress to GPU memory.

[0149] S45, using Segalman formula to calculate Mises stress response power spectral density of element node, and then writing the result to GPU memory, the formula is as follows:

[0150] ;

[0151] Where, is the Mises stress response power spectral density, is the transfer function of the six stress components.

[0152] ;

[0153] is the input power spectral density function.

[0154] S46, calculating node average Mises stress power spectral density by CUDA parallel calculation, each thread accumulating Mises stress power spectral density of associated element node according to corresponding node information, using shared memory to store the accumulation result, and the main thread using arithmetic average method to calculate the average value to get the node average Mises stress, the formula is as follows:

[0155] Assuming that node is shared by 3 elements, then the smoothed stress at node = (stress of element ① at node + stress of element ② at node Stress + element ③ at node (stress) / 3;

[0156] Use thread synchronization (__syncthreads) to ensure consistency of shared memory operations;

[0157] S47. Copy the GPU calculation results from GPU memory back to host memory, write them to the output array, release all GPU memory used by the input and output arrays, and compare the Mises stress response power spectral density curve calculated by the finite element model of the bowl-shaped metal rubber with that of commercial software.

[0158] S5. Using GPU parallel computing technology and frequency domain fatigue analysis methods, multiple frequency domain damage models are used to calculate the random vibration fatigue damage and logarithmic life of unit nodes and global nodes.

[0159] Specifically, it includes the following steps:

[0160] S51. Calculate the root mean square (RMS) stress value at each node using the following formula:

[0161] ;

[0162] in, For root mean square stress, For the first Stress values ​​at each sampling point This represents the number of sampling points.

[0163] Identify the node corresponding to the maximum root mean square stress and regard it as the node most prone to fatigue failure in the structure.

[0164] S52. Flatten the Mises stress response power spectral density of the element node into a continuous matrix, with each column corresponding to an element node and the rows being interpolation frequency points, which facilitates subsequent frequency interpolation and fatigue damage calculation.

[0165] Logarithmically interpolate the Mises stress response power spectral density for each element node to match the frequency points of the input acceleration power spectral density;

[0166] S53. Spectral moments are calculated using numerical integration, as shown in the following formula:

[0167] ;

[0168] in, For the first Spectral moments, For frequency, Power spectral density;

[0169] The zero-point crossing rate is calculated using the following formula:

[0170] ;

[0171] wherein, is the zero crossing rate, is the 2nd order spectral moment, is the 0th order spectral moment.

[0172] The peak crossing rate calculation is performed, with the formula as follows:

[0173] ;

[0174] wherein, is the peak crossing rate, is the 2nd order spectral moment, is the 4th order spectral moment;

[0175] The irregularity factor calculation is performed, with the formula as follows:

[0176] ;

[0177] wherein, is the irregularity factor, is the 0th order spectral moment, is the 2nd order spectral moment, is the 4th order spectral moment;

[0178] The spectral width coefficient calculation is performed, with the formula as follows:

[0179] ;

[0180] wherein, is the peak crossing rate, is the 0th order spectral moment, is the 2nd order spectral moment, is the 4th order spectral moment;

[0181] The above parameters are provided as inputs for the frequency domain fatigue model;

[0182] S54, the fatigue damage of the unit node is calculated by a plurality of frequency domain fatigue models (including Bendat method, Wirsching-Light method, Dirlik method, Tovo-Benasciutti method, Alpha=0.75 method and Zhao-Bake method), and the damage value is stored in the corresponding matrix, and the damage factor is defined as follows:

[0183] ;

[0184] wherein, is the total number of cycles of the constant stress amplitude load block, is the constant stress amplitude of the constant stress amplitude load block, is the constant stress amplitude of the constant stress amplitude load block, The number of cycles to failure under a constant stress amplitude load block, The total number of constant amplitude load blocks when the damage factor is equal to 1 is the failure, and the life calculation criterion is the Basquin formula.

[0185] The specific frequency domain fatigue model formula is as follows:

[0186] Bendat method: ;

[0187] Wirsching-Light method: wherein , , is the slope of the S-N curve.

[0188] Dirlik method: wherein , , , , , , ;

[0189] Tovo-Benasciutti method: wherein is a weight coefficient, and the formula is: ;

[0190] Empirical =0.75 method: ;

[0191] Zhao-Baker method: wherein , , ;

[0192] S55, average damage and life of nodes are calculated by GPU parallel computing, GPU memory is allocated for input data (mapping relationship from nodes to elements, node count, and element node damage data, etc.) and output data (average damage of nodes), and flattened input data is copied from host memory to GPU memory;

[0193] CUDA thread blocks and grids are configured, each thread processes a node, and the grid size is dynamically adjusted to cover all nodes;

[0194] The average damage of nodes is calculated by arithmetic mean, assuming that a node is shared by 3 elements, then the damage of the node after smoothing is = .

[0195] The formula for calculating the logarithmic life is:

[0196] ;

[0197] wherein, is the logarithmic life, is the damage factor.

[0198] S56, the node average damage value calculated by the GPU is copied back from the GPU memory to the host memory, all GPU memories used by the input array and the output array are released, and finally the fatigue damage and the life calculation result are output, and the damage result cloud map and the logarithmic life cloud map of the bowl-shaped metal rubber finite element model obtained by calculation are shown in Figure 9 .

[0199] The present application combines random vibration theory and finite element method, simultaneously introduces MPI, OpenMP and GPU parallel technology to realize high-performance calculation, so as to solve the related problems of metal rubber random vibration fatigue, and can quickly realize the simulation calculation of large-scale random vibration fatigue damage of metal rubber, and the heterogeneous parallel acceleration ratio, and the time consumption comparison of serial and parallel calculation under the same model.

Claims

1. A heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber, characterized in that, Includes the following steps: S1. Input the finite element model data of large-scale metal-rubber components; S2. Using CPU-based multi-threaded parallel technology, calculate the total stiffness matrix and total mass matrix of the finite element model; S3. Combining MPI multi-process parallel technology and OpenMP multi-threaded parallel technology, calculate the frequency response function of the finite element model; S4. Using GPU parallel computing technology, calculate random vibration response based on frequency response function, including unit node Mises stress power spectral density and node average Mises stress power spectral density. S5. Using GPU parallel computing technology and frequency domain fatigue analysis methods, multiple frequency domain damage models are used to calculate the random vibration fatigue damage and logarithmic life of unit nodes and global nodes.

2. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 1, characterized in that, S2 specifically includes: S21. Each CPU thread is allocated an independent memory space through a multi-threaded parallel framework, and the elements in the finite element mesh are allocated to multiple threads for parallel computation through dynamic scheduling. S22. Within the same parallel region, each thread computes the element stiffness matrix and element mass matrix of the assigned element in parallel using OpenMP multithreading, based on the assigned element. S23. Merge the local Triplet lists of each thread into the global Triplet list through the critical section, and assemble the element matrix into the corresponding positions of the total stiffness matrix and the total mass matrix through it. S24. Apply boundary conditions using the large number method and update the total stiffness matrix and total mass matrix.

3. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 2, characterized in that, In step S22, the calculation method for the element stiffness matrix after B-Bar correction is as follows: ； The element mass matrix is ​​in lumped mass matrix form and is calculated using the following formula: ； in, The element stiffness matrix, The unit mass matrix, The geometric matrix after B-Bar correction. For material matrix, It is a unit shape function matrix.

4. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 1, characterized in that, S3 specifically includes: S31. Initialize the MPI distributed computing environment. The main process distributes the total stiffness matrix and total mass matrix to the other processes through MPI multi-process broadcast. S32. Call the MUMPS open-source distributed solver to perform LU decomposition on the total stiffness matrix; S33. Call the PARPACK open-source library to collaboratively calculate the system's eigenvalues ​​and eigenvectors using the shift-inversion method; S34. Calculate the displacement response power spectral density of the nodes.

5. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 4, characterized in that, The shift-reversal method is used to solve the generalized eigenvalue problem. Transform into ,in , To find the smallest eigenvalue Transform into the largest eigenvalue Perform a fast solution.

6. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 1, characterized in that, S4 specifically includes: S41. Preprocess the input data and copy it from host memory to GPU memory; S42. Configure CUDA thread blocks and mesh, with each thread processing one cell, local node, and frequency point; S43. Calculate the strain using the B-Bar corrected B matrix, then calculate the Gaussian integral point stress using the material D matrix, and finally calculate the element nodal stress using the shape function matrix. S44. Calculate the Mises stress response power spectral density of the element node using the Segalman formula; S45. Calculate the node-average Mises stress power spectral density in parallel using CUDA; S46. Copy the GPU computation results from GPU memory back to host memory.

7. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 1, characterized in that, S5 specifically includes: S51. Spectral moments are calculated through numerical integration, and zero-point crossing rate, peak crossing rate, and irregularity factor are further calculated. S52. Calculate the fatigue damage of unit nodes using multiple frequency domain fatigue models; S53. The average damage and lifetime of nodes are calculated in parallel using GPU, where the average damage of nodes is calculated using an arithmetic mean. S54. Output fatigue damage and life calculation results.

8. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 7, characterized in that, The various frequency domain fatigue models include the Bendat method, the Wirsching-Light method, the Dirlik method, the Tovo-Benasciutti method, the Alpha=0.75 method, and the Zhao-Bake method.

9. The heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 1, characterized in that, In S4 and S5, both methods utilize GPU parallel computing technology and include the steps of flattening data such as node displacement frequency response function and B-Bar matrix into continuous arrays, copying them from host memory to GPU memory, and configuring CUDA thread blocks and meshes.

10. A heterogeneous parallel method for predicting random vibration damage and fatigue life of metal rubber according to claim 7, characterized in that, The formula for calculating the logarithmic lifetime is: ； in, For logarithmic lifetime, It is a damage factor.

Images