Multi-layer composite material damage analysis method based on numerical simulation calculation

Abstract

The invention discloses a multilayer composite material damage analysis method based on numerical simulation calculation, and relates to the technical field of composite material damage analysis, and the method comprises the following steps: S1, constructing an initial single-layer stiffness matrix C'of each layer; s2, constructing a single-layer stiffness matrix Alayer and a single-layer mass matrix Mlayer of each layer; s3, setting a free boundary condition; s4, constructing a global stiffness matrix (Aglobal) and a global quality matrix (Mglobal); s5, the global stiffness matrix Aglobal is subjected to regularization processing, and # imgabs0 # is obtained; and S6, constructing a frequency dispersion equation, solving the frequency dispersion equation to obtain a characteristic value lambda and a characteristic vector phi, and calculating a phase velocity vp and a group velocity vg of the Lamb wave. The method provided by the invention is used for the damage analysis of the orthotropic multilayer composite material, and the calculation accuracy is improved.

Description

Technical Field

[0001] The present invention relates to the technical field of composite material damage analysis, and in particular to a multilayer composite material damage analysis method based on numerical simulation calculation. Background Art

[0002] As a nondestructive testing method, Lamb wave is widely used in the field of health monitoring of composite materials due to its sensitivity to internal damage of materials, such as nondestructive monitoring through Lamb wave and damage analysis of composite materials through numerical simulation calculations.

[0003] Composite material damage analysis through numerical simulation calculation means first modeling the composite material, then generating a wave number k through numerical simulation (rather than actually applying a Lamb wave to the composite material), and then calculating and generating the propagation characteristic parameters of the Lamb wave in the composite material on this basis for damage analysis of the composite material.

[0004] At present, the damage analysis of composite materials calculated by numerical simulation is mainly focused on single-layer materials. However, for multi-layer composite materials, the interlayer coupling effect is often ignored. In multi-layer composite materials, the interlayer coupling has a significant impact on the propagation of Lamb waves. Ignoring this may lead to errors.

[0005] In view of this, this invention is proposed. Summary of the invention

[0006] The purpose of the present invention is to solve the shortcomings existing in the prior art and to propose a multilayer composite material damage analysis method based on numerical simulation calculation, which is used for damage analysis of orthotropic multilayer composite materials. The propagation characteristic parameters of Lamb waves in the multilayer composite materials are calculated by numerical simulation calculation method, and the accuracy of calculation is improved.

[0007] In order to achieve the above object, the present invention adopts the following technical solutions: A multilayer composite material damage analysis method based on numerical simulation calculation includes the following steps: S1: Construct the initial single-layer stiffness matrix C' of each layer; S2: Construct the single-layer stiffness matrix A of each layer according to the initial single-layer stiffness matrix C' of each layer layer and the single-layer mass matrix M layer ; S3: Set free boundary conditions; S4: Accumulate layer by layer to obtain the global stiffness matrix A global and the global mass matrix M global , in the process of layer-by-layer accumulation and construction, inter-layer continuous processing and calculation are performed according to the inter-layer coupling effect; S5: For the global stiffness matrix A global After regularization, we can obtain ; S6: Based on and M global Construct the dispersion equation, solve the dispersion equation to obtain the eigenvalue λ and eigenvector Φ, and calculate the group velocity v of the Lamb wave p and phase velocity v g .

[0008] Further, the S1 comprises the following steps: S11: Construct the original single-layer stiffness matrix C in the following form: , in, , , , , , , , , , in, , , , are the Poisson's ratios of the xy, yx, xz, and yz planes, respectively. , , are the elastic moduli in the x, y, and z directions, respectively, , , denote the shear modulus in the xy, yz, and xz planes, respectively; S12: Construct the rotation matrix T in the following form: , in, , , Where θ is the rotation angle; S13: Calculate the initial single layer stiffness matrix C' according to the following formula .

[0009] Further, the S2 comprises the following steps: S21: Using the Gauss-Legendre integral method, the single-layer stiffness matrix A of each layer is calculated by the following formula layer and the single-layer mass matrix M layer The elements in: , , Where N is the number of integration points, i, j = 0, 1, ..., N-1, a and b are the index values ​​of the Gauss-Legendre integration points, is the N-1 order Legendre polynomial basis function, is the weight, h is the layer thickness, is the density.

[0010] Further, the S3 comprises the following steps: S31: For the top layer of a multilayer composite, the free boundary condition is set by the following formula: , , S32: For the bottom layer of a multilayer composite material, the free boundary condition is set by the following formula: , , Wherein, M is the number of layers of the multilayer composite material.

[0011] Further, the S4 comprises the following steps: S41: The global stiffness matrix of the first layer is accumulated according to the following formula: and the global mass matrix , , , S42: The global stiffness matrix of the next layer is accumulated according to the following formula: and the global mass matrix , , , Where, c=1,2,...M, M is the number of layers of the multilayer composite material; S43: By superimposing adjacent layer matrix blocks A upper and A lower , update the global stiffness matrix and the global mass matrix The elements satisfy the displacement and stress continuity conditions; S44: Execute S42 and S43 repeatedly until the cumulative construction is completed and the global stiffness matrix A is obtained. global and the global mass matrix M global .

[0012] Further, the S43 comprises the following steps: S431: Construct the matrix block A at the bottom of the previous layer upper , and calculate A according to the following formula upper Elements in: , Where F is the number of basis functions used in discretization; S432: Build the matrix block A at the top of the current layer lower , and calculate A according to the following formula lower Elements in: , S433: Update according to the following formula Elements in: .

[0013] Further, the S5 comprises the following steps: S51: calculated by the following formula , , in, =10 -6 , I is in A global The identity matrix of the same order, A global condition number.

[0014] Further, the S6 comprises the following steps: S61: Construct the dispersion equation in the following form: , S62: Solve the dispersion equation and obtain a series of eigenvalues ​​λ= , and the eigenvector Ф corresponding to each one, where ω is the angular frequency; S63: Delete the imaginary eigenvalue λ, and obtain the effective eigenvalue λ and the effective eigenvector Φ; S64: Calculate the phase velocity v of the Lamb wave according to the following formula p and phase velocity v g， , , Where k is the wave number of the multilayer composite material and f is the frequency.

[0015] Further, the method further comprises the following steps: S7: With frequency f as x-axis and phase velocity v as p Or phase velocity v g As the y-axis, a visualization image is constructed with the effective eigenvector Ф as the z-axis.

[0016] Compared with the prior art, the present invention has the following beneficial effects: It can realize the calculation of the phase velocity and group velocity of Lamb waves in different frequency ranges of orthogonal anisotropic multilayer composite materials; by accurately modeling the global stiffness matrix and global mass matrix of multilayer composite materials, and performing interlayer continuous processing calculations based on the interlayer coupling effect, it can accurately reflect the stress-strain relationship of each layer. Compared with the existing technology, the accuracy of Lamb wave propagation characteristic analysis is improved, making the damage analysis and research of multilayer composite materials more accurate and reliable, and can provide accurate theoretical support for damage detection of multilayer composite materials. BRIEF DESCRIPTION OF THE DRAWINGS

[0017] Figure 1 It is a flowchart of a damage analysis method for multilayer composite materials based on numerical simulation calculation; Figure 2 A schematic diagram of the visualization image. DETAILED DESCRIPTION

[0018] The technical solutions in the embodiments of the present invention will be described clearly and completely below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all the embodiments.

[0019] Embodiment 1: A damage analysis method for multilayer composite materials based on numerical simulation calculations, such as Figure 1 As shown, the following steps are included: S1: Construct the initial single-layer stiffness matrix C' of each layer, S2: Construct the single-layer stiffness matrix A of each layer according to the initial single-layer stiffness matrix C' of each layer layer and the single-layer mass matrix M layer， S3: Set free boundary conditions, S4: Accumulate layer by layer to obtain the global stiffness matrix A global and the global mass matrix M global In the process of layer-by-layer accumulation and construction, continuous inter-layer processing and calculation are performed according to the inter-layer coupling effect. S5: For the global stiffness matrix A global After regularization, we can obtain , S6: Based on and M global Construct the dispersion equation, solve the dispersion equation to obtain the eigenvalue λ and eigenvector Φ, and calculate the group velocity v of the Lamb wave p and phase velocity v g .

[0020] The multilayer composite material damage analysis method based on numerical simulation calculation in this embodiment can realize the calculation of the phase velocity and group velocity of Lamb waves in different frequency ranges of orthogonal anisotropic multilayer composite materials; by accurately modeling the global stiffness matrix and global mass matrix of the multilayer composite material, and performing interlayer continuous processing calculation according to the interlayer coupling effect, it can accurately reflect the stress-strain relationship of each layer, and improve the accuracy of Lamb wave propagation characteristic analysis compared with the prior art, making the damage analysis and research of multilayer composite materials more accurate and reliable, and can provide accurate theoretical support for damage detection of multilayer composite materials.

[0021] In an optional embodiment, the S1 comprises the following steps: S11: Construct the original single-layer stiffness matrix C in the following form: , in, , , , , , , , , , in, , , , are the Poisson's ratios of the xy, yx, xz, and yz planes, respectively. , , They are the elastic modulus in the x, y, and z directions, respectively, which reflect the material's ability to resist elastic deformation in the x, y, and z directions. , , They represent the shear modulus of the xy, yz, and xz planes, respectively, which reflect the material's ability to resist shear deformation in the xy, yz, and xz planes. In addition, Poisson's ratio, elastic modulus, and shear modulus are material properties, where It can be calculated by the following formula: .

[0022] S12: Construct the rotation matrix T in the following form: , in, , , Where θ is the rotation angle, which is used to describe the degree of rotation of the material layer relative to the reference coordinate system.

[0023] S13: Calculate the initial single-layer stiffness matrix C' according to the following formula: .

[0024] In this optional embodiment, the original single-layer stiffness matrix C is constructed based on the orthogonal anisotropy of the material, and the initial single-layer stiffness matrix C' is constructed by rotation transformation based on the rotation angle θ of each layer, thereby ensuring that the matrix C' can accurately reflect the stress-strain relationship of each layer.

[0025] In an optional embodiment, the S2 comprises the following steps: S21: Using the Gauss-Legendre integral method, the single-layer stiffness matrix A of each layer is calculated by the following formula layer and the single-layer mass matrix M layer The elements in: , , Where N is the number of integration points, i, j = 0, 1, ..., N-1, a and b are the index values ​​of the Gauss-Legendre integration points, is the N-1 order Legendre polynomial basis function, is the weight, h is the layer thickness, is the density; additionally, density is a material property.

[0026] Preferably, according to the convergence analysis, the number of integration points N can be 5, at which time the relative error is less than 1%.

[0027] In an optional embodiment, S3 includes the following steps: S31: For the top layer of a multilayer composite, the free boundary condition is set by the following formula: , .

[0028] S32: For the bottom layer of a multilayer composite material, the free boundary condition is set by the following formula: , , Wherein, M is the number of layers of the multilayer composite material.

[0029] In an optional embodiment, the S4 comprises the following steps: S41: The global stiffness matrix of the first layer is accumulated according to the following formula: and the global mass matrix , , .

[0030] S42: The global stiffness matrix of the next layer is accumulated according to the following formula: and the global mass matrix , , , Wherein, c=1,2,...M, and M is the number of layers of the multilayer composite material.

[0031] S43: By superimposing adjacent layer matrix blocks A upper and A lower , update the global stiffness matrix and the global mass matrix elements satisfy the displacement and stress continuity conditions.

[0032] S44: Execute S42 and S43 repeatedly until the cumulative construction is completed, that is, the cumulative calculation reaches and , obtain the global stiffness matrix A global and the global mass matrix M global .

[0033] In this optional embodiment, first according to and calculate and , and then accumulate the calculation layer by layer until the cumulative calculation is completed and .

[0034] In an optional embodiment, the S43 includes the following steps: S431: Construct the matrix block A at the bottom of the previous layer upper , and calculate A according to the following formula upper Elements in: , Here, F is the number of basis functions used in discretization.

[0035] S432: Build the matrix block A at the top of the current layer lower , and calculate A according to the following formula lower Elements in: .

[0036] S433: Update according to the following formula Elements in: .

[0037] Starting from the second layer of the multilayer composite material, it is necessary to ensure the stress continuity between two adjacent layers. In this optional embodiment, the stress continuity between two adjacent layers can be ensured by performing interlayer continuity processing calculation according to the interlayer coupling effect.

[0038] Specifically, by diagonally splicing matrix blocks and combining interlayer continuity processing, the global stiffness matrix A is ensured to be global and the global mass matrix M global The accuracy and completeness of the modeling of the overall characteristics of multilayer composite materials are achieved, providing a solid foundation for solving the dispersion equation.

[0039] In an optional embodiment, the S5 comprises the following steps: S51: calculated by the following formula , , in, =10 -6 , I is in A global The identity matrix of the same order, A global In addition, it should be noted that when 10 −8 ≤ ≤10 −4 When the phase velocity is calculated, the deviation is less than 0.5%.

[0040] In an optional embodiment, the S6 comprises the following steps: S61: Construct the dispersion equation in the following form: , S62: Solve the dispersion equation and obtain a series of eigenvalues ​​λ= , and the eigenvector Φ corresponding to each one, where ω is the angular frequency.

[0041] In this optional embodiment, the dispersion equation can be solved by numerical methods (such as eigh function).

[0042] Different angular frequencies ω correspond to different propagation modes of Lamb waves in multilayer composite materials, such as the symmetric mode S n and antisymmetric mode A n ; The eigenvector Ф describes the vibration mode of the Lamb wave at the corresponding angular frequency. The eigenvector Ф contains the displacement components u1, u2 and u3 of the Lamb wave and can be expressed as follows: , The amplitude calculation formula is as follows: .

[0043] S63: Delete the imaginary eigenvalue λ, and obtain the effective eigenvalue λ and the effective eigenvector Φ.

[0044] Since the imaginary eigenvalue λ corresponds to an attenuated wave or a non-physical solution, in this optional embodiment, the imaginary eigenvalue λ is deleted, that is, the real eigenvalue λ is screened out from the series of eigenvalues ​​λ calculated and obtained in S62, and then the effective eigenvalue λ and the effective eigenvector Φ are obtained to characterize the propagating Lamb wave mode.

[0045] S64: Calculate the phase velocity v of the Lamb wave according to the following formula p and phase velocity v g， , , Where k is the wave number of the multilayer composite material and f is the frequency.

[0046] For the wave number k, an initial wave number guess value k0 can be given first. This initial value can be selected based on prior knowledge or experience of the physical problem, and then the Newton-Raphson method is used to iteratively solve and obtain the wave number k. In addition, the angular frequency ω and the frequency f can be converted according to the following formula.

[0047] In this optional embodiment, by solving the generalized eigenvalue problem, combined with eigenvalue deletion and screening, the dispersion equation can be accurately solved and the phase velocity v of the Lamb wave can be calculated. p and phase velocity v g , and then comprehensively describe the propagation characteristics of Lamb waves in multilayer composite materials.

[0048] In an optional embodiment, the method further includes the following steps: S7: With frequency f as x-axis and phase velocity v as p Or phase velocity v gThe visual image is constructed with the amplitude of the effective eigenvector Ф as the y-axis and the amplitude of the effective eigenvector Ф as the z-axis.

[0049] In this optional embodiment, the frequency f is the x-axis and the phase velocity v is the p Or phase velocity v g The visual image is constructed with the amplitude of the effective eigenvector Ф as the y-axis and the amplitude of the effective eigenvector Ф as the z-axis, which can clearly show the changing trends of the phase velocity and group velocity with the change of frequency, making the analysis results of the propagation characteristics of Lamb waves more intuitive and easy to understand, which is convenient for technical personnel to observe and analyze, and improves the efficiency and accuracy of damage analysis and research.

[0050] Specifically, phase velocity refers to the propagation speed of the phase front of the lamb wave, so the curve constructed with phase velocity can show the phase propagation speed of the wave at different frequencies; group velocity refers to the energy propagation speed of the lamb wave, so the dispersion curve constructed with group velocity can show the energy propagation speed of the wave at different frequencies; by observing the dispersion curve, the Lamb wave mode at different frequencies can be determined; the dispersion phenomenon can be identified through the bending and changes of the dispersion curve, which is of great significance in non-destructive testing and structural health testing.

[0051] Specifically, To further demonstrate the calculation results of this embodiment, a carbon fiber reinforced material (CFRP) T700M913 is used, with a size of 500*500mm, 8 layers, and the rotation angle of the upper four layers is [0 / 90] 2s , the rotation angle of the next four layers is [90 / 0] 2s , the thickness of each layer is 0.15 mm as an example for demonstration calculation. The material properties of the material are shown in Table 1.

[0052] Table 1 Material properties of carbon fiber reinforced materials <![CDATA[ρ(kg / m 3 )]]> <![CDATA[ E x ( GPα )]]> <![CDATA[ E y ( GPα )]]> <![CDATA[ E z ( GPα )]]> <![CDATA[ G xy ( GPα )]]> 1.56 135 10 10 5 <![CDATA[ G xz ( GPα )]]> <![CDATA[ G yz ( GPα )]]> <![CDATA[ v xy ]]> <![CDATA[ v xz ]]> <![CDATA[ v yz ]]> 5 5 0.3 0.25 0.25

[0053] Visualization images such as Figure 2 As shown, the color bar on the right represents the amplitude of the effective eigenvector Φ, and this embodiment is demonstrated by taking the ordinate as the phase velocity as an example.

[0054] The above description is only a preferred specific implementation manner of the present invention, but the protection scope of the present invention is not limited thereto. Any technician familiar with the technical field can make equivalent replacements or changes according to the technical scheme and inventive concept of the present invention within the technical scope disclosed by the present invention, which should be covered by the protection scope of the present invention.

Claims

1. A damage analysis method for multilayer composite materials based on numerical simulation calculation, characterized in that: The steps include: S1: Construct the initial single-layer stiffness matrix C' of each layer; S2: Construct the single-layer stiffness matrix A of each layer according to the initial single-layer stiffness matrix C' of each layer layer and the single-layer mass matrix M layer ; S3: Set free boundary conditions; S4: Accumulate layer by layer to obtain the global stiffness matrix A global and the global mass matrix M global , in the process of layer-by-layer accumulation and construction, inter-layer continuous processing and calculation are performed according to the inter-layer coupling effect; S5: For the global stiffness matrix A global After regularization, we can obtain ; S6: Based on and M global Construct the dispersion equation, solve the dispersion equation to obtain the eigenvalue λ and eigenvector Ф, and calculate the phase velocity v of the Lamb wave p and group velocity v g .

2. A multilayer composite material damage analysis method based on numerical simulation calculation according to claim 1, characterized in that: The S1 comprises the following steps: S11: Construct the original single-layer stiffness matrix C in the following form: ， in, ， ， ， ， ， ， ， ， ， in, , , , are the Poisson's ratios of the xy, yx, xz, and yz planes, respectively. , , are the elastic moduli in the x, y, and z directions, respectively, , , denote the shear modulus in the xy, yz, and xz planes, respectively; S12: Construct the rotation matrix T in the following form: ， in, ， ， Where θ is the rotation angle; S13: Calculate the initial single-layer stiffness matrix C' according to the following formula: 。 3. A multilayer composite material damage analysis method based on numerical simulation calculation according to claim 1, characterized in that: The S2 comprises the following steps: S21: Using the Gauss-Legendre integral method, the single-layer stiffness matrix A of each layer is calculated by the following formula layer and the single-layer mass matrix M layer The elements in: ， ， Where N is the number of integration points, i, j = 0, 1, ..., N-1, a and b are the index values ​​of the Gauss-Legendre integration points, is the N-1 order Legendre polynomial basis function, is the weight, h is the layer thickness, is the density.

4. A multilayer composite material damage analysis method based on numerical simulation calculation according to claim 3, characterized in that: The S3 comprises the following steps: S31: For the top layer of a multilayer composite, the free boundary condition is set by the following formula: ， ， S32: For the bottom layer of a multilayer composite material, the free boundary condition is set by the following formula: ， ， Wherein, M is the number of layers of the multilayer composite material.

5. The damage analysis method of multilayer composite materials based on numerical simulation calculation according to claim 1, characterized in that: The S4 comprises the following steps: S41: The global stiffness matrix of the first layer is accumulated according to the following formula: and the global mass matrix , ， ， S42: The global stiffness matrix of the next layer is accumulated according to the following formula: and the global mass matrix , ， ， Where, c=1,2,...M, M is the number of layers of the multilayer composite material; S43: By superimposing adjacent layer matrix blocks A upper and A lower , update the global stiffness matrix and the global mass matrix The elements satisfy the displacement and stress continuity conditions; S44: Execute S42 and S43 repeatedly until the cumulative construction is completed and the global stiffness matrix A is obtained. global and the global mass matrix M global .

6. A multilayer composite material damage analysis method based on numerical simulation calculation according to claim 5, characterized in that: The S43 comprises the following steps: S431: Construct the matrix block A at the bottom of the previous layer upper , and calculate A according to the following formula upper Elements in: ， Where F is the number of basis functions used in discretization; S432: Build the matrix block A at the top of the current layer lower , and calculate A according to the following formula lower Elements in: ， S433: Update according to the following formula Elements in: 。 7. The damage analysis method of multilayer composite materials based on numerical simulation calculation according to claim 1, characterized in that: The S5 comprises the following steps: S51: calculated by the following formula , ， in, =10 -6 , I is in A global The identity matrix of the same order, A global condition number.

8. The damage analysis method of multilayer composite materials based on numerical simulation calculation according to claim 1, characterized in that: The S6 comprises the following steps: S61: Construct the dispersion equation in the following form: ， S62: Solve the dispersion equation and obtain a series of eigenvalues ​​λ= , and the eigenvector Ф corresponding to each one, where ω is the angular frequency; S63: Delete the imaginary eigenvalue λ, and obtain the effective eigenvalue λ and the effective eigenvector Φ; S64: Calculate the phase velocity v of the Lamb wave according to the following formula p and phase velocity v g， ， ， Where k is the wave number of the multilayer composite material and f is the frequency.

9. A multilayer composite material damage analysis method based on numerical simulation calculation according to claim 8, characterized in that: The following steps are also included: S7: With frequency f as x-axis and phase velocity v as p Or phase velocity v g As the y-axis, a visualization image is constructed with the effective eigenvector Ф as the z-axis.

Images