An ultrasonic guided wave dispersion curve calculation method

Abstract

The application provides an ultrasonic guided wave dispersion curve calculation method, and belongs to the technical field of ultrasonic measurement. With the increase of the thickness of a to-be-measured plate structure and the increase of the complexity of the structure, if a traditional finite element analysis method is used when a finite element model of the to-be-measured plate structure is constructed, a high-order polynomial interpolation function needs to be used, and a large number of units need to be divided, so that the number of finite element grids increases, and the calculation efficiency is reduced. In the application, the shape function of an interval B-spline small element is used as an interpolation function for solving, the calculation accuracy can be ensured, a large number of units do not need to be divided, the number of finite element grids is reduced, and the calculation efficiency can be improved.

Description

TECHNICAL FIELD

[0001] The application relates to an ultrasonic guided wave dispersion curve calculation method and belongs to the technical field of ultrasonic wave measurement. BACKGROUND

[0002] Ultrasonic guided wave detection is a new nondestructive testing technology developed from traditional ultrasonic detection technology. Compared with traditional ultrasonic bulk waves, guided waves only need single-point excitation to complete the detection of the entire measured medium, and guided waves have the advantages of fast propagation speed, small energy attenuation and high detection efficiency, and are widely used in rapid detection of plate structures. The propagation characteristics of guided waves in different waveguide media contain some basic information of ultrasonic guided waves, and the phase velocity and group velocity dispersion curves of ultrasonic wave propagation are mainly solved by establishing a related model, thereby providing theoretical support for the selection of guided wave frequency and the selection of probe oblique incidence angle in actual detection.

[0003] At present, the dispersion curve solving method of a plate structure is mainly to obtain the relationship between frequency and wave number by solving a characteristic equation and a finite element method, and then draw the phase velocity and group velocity dispersion curves of the measured plate structure. When the finite element method is used, the finite element grid is usually divided, and a high-order polynomial basis function is used for solving. Since the size of the finite element grid is generally controlled at λ / 20~λ / 10 (λ is the minimum wavelength), with the increase of the thickness of the plate structure, the number of finite element grids increases, resulting in low calculation efficiency of the solving process. SUMMARY

[0004] The application aims to provide an ultrasonic guided wave dispersion curve calculation method, which is used to solve the problem of low calculation efficiency of solving dispersion curves caused by the increase of plate thickness.

[0005] In order to achieve the above-mentioned purpose, the application provides an ultrasonic guided wave dispersion curve calculation method, which comprises the following steps:

[0006] S1, establishing a finite element model in the thickness direction of a plate structure to be measured, and dividing the finite element model into units;

[0007] S2, calculating the unit displacement of each unit in the finite element model by using the shape function constructed by the interval B-spline small element;

[0008] S3, calculating the unit strain of the corresponding unit according to the unit displacement calculated in step S2;

[0009] S4, calculating the stiffness matrix and the mass matrix of the corresponding unit according to the unit displacement and the unit strain of each unit;

[0010] S5, according to the stiffness matrix and the mass matrix of each unit, a wave equation is established, and the wave equation is solved, the relationship between the wave number and the frequency of the ultrasonic guided wave is obtained, and then the dispersion curve of the ultrasonic guided wave when transmitting in the plate structure to be measured is drawn.

[0011] With the increase of the thickness of the plate structure to be measured, when constructing the finite element model of the plate structure to be measured, if the traditional finite element analysis method is used, a high-order polynomial interpolation function needs to be used, and a large number of units need to be divided, so that the number of finite element grids increases, and the calculation efficiency is reduced. The shape function of the interval B-spline wavelet element is used as the interpolation function for solving, which can ensure the calculation accuracy, and does not need to divide a large number of units, the number of finite element grids is reduced, and the calculation efficiency can be improved.

[0012] Further, in the above method, the shape function of the interval B-spline wavelet element constructed in step S2 is obtained by the following formula:

[0013]

[0014] In the formula, u is the shape function of the interval B-spline wavelet element, Φ is the row vector of the scale function of the interval B-spline wavelet function, T e is the element conversion matrix of the one-dimensional interval B-spline wavelet element.

[0015] A specific method is provided to obtain the shape function of the interval B-spline wavelet element, which is convenient to realize.

[0016] Further, in the above method, in step S2, the element displacement of the e-th element in the finite element is calculated by the following formula:

[0017]

[0018] In the formula, u (e) is the element displacement, is the shape function of the interval B-spline wavelet element, q (e) is the node displacement, e i (κx-ωt) is a simple harmonic wave, i is an imaginary symbol, κ is the wave number of the ultrasonic guided wave, x represents the propagation direction of the ultrasonic guided wave, ω is the frequency of the ultrasonic guided wave, and t is time.

[0019] A specific method is provided to calculate the element displacement, which is convenient to realize.

[0020] Further, in the above method, in step S3, the element strain of each unit is calculated by the following formula:

[0021]

[0022]

[0023] In the formula, ε e is the unit strain, u (e) is the unit displacement.

[0024] A specific method is provided to calculate the unit strain, facilitating implementation.

[0025] Further, in the above method, in step S4, the stiffness matrix and the mass matrix of each unit are obtained by substituting the unit displacement and the unit strain into the discrete form of the Hamilton formula; the discrete form of the Hamilton formula is:

[0026]

[0027] In the formula, t0 and t1 respectively represent the starting time and the ending time of the motion of the mass point, n el is the number of units of the finite element, e = 1, 2,..., n el , ε (e) represents the unit strain, δ(ε (e)T ) represents the virtual strain, δ(u (e)T ) represents the virtual displacement, C e represents the material stiffness matrix, u (e) represents the unit displacement, ρ e represents the material density of the plate structure to be measured, represents the second-order inverse of the unit displacement with respect to time t.

[0028] A specific calculation method is provided to solve the stiffness matrix and the mass matrix, facilitating implementation of the present application.

[0029] Further, in the above method, after substituting the unit displacement and the unit strain into the discrete form of the Hamilton formula, the following formula is obtained, and the stiffness matrix and the mass matrix of the corresponding unit are obtained by solving the formula;

[0030]

[0031] In the formula, δq (e)T represents the virtual displacement, and are three stiffness matrices in the corresponding unit, m (e) is the mass matrix of the corresponding unit, κ represents the wave number of the ultrasonic guided wave, ω represents the frequency of the ultrasonic guided wave, q (e) represents the node displacement, and i is the imaginary symbol.

[0032] Further, in the above method, in step S5, when establishing the wave equation, the specific operator symbols Γ 0,0 , Γ 1,0 , Γ1,0 and Γ 1,1 simplifying the stiffness matrix and the mass matrix;

[0033]

[0034]

[0035] wherein T e represents a unit conversion matrix, l e represents the length of the e-th unit, and Φ represents a row vector of a scaling function of an interval B-spline wavelet function.

[0036] The stiffness matrix and the mass matrix are simplified by introducing specific operators, thereby reducing the amount of calculation in the subsequent calculation process and improving the calculation efficiency.

[0037] Further, in the above method, the stiffness matrix and the mass matrix are integrated into a global stiffness matrix and a global mass matrix by the following formula, and then the wave equation is established by the global stiffness matrix and the global mass matrix.

[0038]

[0039] wherein K1, K2 and K3 are global stiffness matrices corresponding to and M is a global mass matrix corresponding to m (e .

[0040] Further, in the above method, the wave equation is:

[0041] [K1+iκK2+κ 2 K3-ω 2 M] O U=0

[0042] wherein κ represents the wave number of the ultrasonic guided wave, i represents the imaginary symbol, K1, K2 and K3 are global stiffness matrices corresponding to and M is a global mass matrix corresponding to m (e) , and U represents the displacement vector at the node, and subscript O is the total number of degrees of freedom.

[0043] Further, in the above method, when solving the wave equation, a diagonal matrix T0 is introduced to eliminate the imaginary symbol i in the wave equation, and the diagonal matrix T0 is:

[0044]

[0045] wherein subscript O is the total number of degrees of freedom.

[0046] When solving the wave equation, introducing a diagonal matrix T0 to eliminate the imaginary part i in the wave equation can reduce the amount of computation and improve computational efficiency. Attached Figure Description

[0047] Figure 1 This is a flowchart illustrating the method for calculating the dispersion curve of ultrasonic guided waves in an embodiment of the present invention.

[0048] Figure 2 This is a schematic diagram of the 4th order 3-scale interval B-spline wavelet scaling function in an embodiment of the method of the present invention;

[0049] Figure 3 This is a schematic diagram showing the correspondence between the number of nodes and node coordinates within a unit in an embodiment of the method of the present invention;

[0050] Figure 4 This is a semi-analytical wavelet finite element model of the plate structure in the embodiment of the invention.

[0051] Figure 5 This is a comparison diagram of the dispersion curves of anisotropic homogeneous plate structures obtained by using the present invention and the conventional finite element method in the embodiments of the present invention.

[0052] Figure 6 This is a time comparison diagram of solving the dispersion curve of an isotropic mean plate structure using the present invention and using the conventional finite element method in the embodiments of the present invention;

[0053] Figure 7 This is a schematic diagram of a composite material laminate model in an embodiment of the method of the present invention;

[0054] Figure 8 This is a comparison diagram of the phase velocity dispersion curves of composite laminates solved using the present invention and the conventional finite element method in the embodiments of the present invention.

[0055] Figure 9 This is a comparison diagram of the velocity dispersion curves of composite laminate groups solved using the present invention and the conventional finite element method in the embodiments of the present invention.

[0056] Figure 10 This is a time comparison graph showing the solution of the dispersion curve of composite laminate using the present invention and the conventional finite element method in the embodiments of the present invention. Detailed Implementation

[0057] To make the objectives, technical solutions, and advantages of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments.

[0058] Method Implementation Examples:

[0059] like Figure 1As shown, the ultrasonic guided wave dispersion curve calculation method of the present application comprises the following steps:

[0060] 1. A semi-analytical method is used to establish a finite element model of the plate structure to be measured, and the finite element model has the same size as the plate structure to be measured in the thickness direction. At the same time, according to the requirement of solution accuracy, different scale interval B-spline wavelet functions are selected, and the scale function of the interval B-spline wavelet element and the element conversion matrix of the one-dimensional interval B-spline wavelet element are obtained. The shape function of the one-dimensional interval B-spline wavelet element is constructed by using the scale function of the interval B-spline wavelet element and the element conversion matrix of the one-dimensional interval B-spline wavelet element.

[0061] As shown in Figure 2 is a semi-analytical wavelet finite element model of the plate structure to be measured, X direction is the wave propagation direction of the ultrasonic guided wave, y-z is the cross section of the plate structure to be measured, and y direction is infinite wide. The wave number of the ultrasonic guided wave is κ, and the frequency is ω.

[0062] In this embodiment, according to the requirement of solution accuracy, a 4th order 3-scale interval B-spline wavelet function is selected. The scale function of the 4th order 3-scale interval B-spline wavelet function has 11 functions, which are respectively: 0 boundary scale function 1 boundary scale function and internal scale function and

[0063] The interval B-spline wavelet (BSWI) function is constructed by giving a simple node sequence and connecting them together using a piecewise polynomial function between the nodes, so as to obtain a certain order of overall smoothness. For any one-dimensional function f(x) in the solution interval [a, b], it can be converted to the standard solution interval [0, 1] by linear mapping relationship ξ=(x-a) / (b-a), so it only needs to construct m order B-spline wavelet function in the solution interval [0, 1]. In order to ensure that there is at least one internal wavelet in the interval [0, 1], the following conditions must be met:

[0064] 2 j ≥2m-1

[0065] Where m and j are the order and scale parameter of the interval B-spline wavelet function, respectively. By using the m order 0 scale interval B-spline wavelet scale function, the interval B-spline wavelet scale function of any j scale can be obtained by the following recursive formula The recursive formula is:

[0066]

[0067] Therefore, the row vector of the scaling function of the interval B-spline wavelet function on the interval [0, 1] can be expressed as:

[0068]

[0069] where the value range of ξ is [0, 1]. As shown in the following formula, it is the scaling function of the 4-order 3-scale interval B-spline wavelet function. Figure 3

[0070] When the m-order j-scale interval B-spline wavelet function is used as an interpolation function construction unit, the number of nodes on the unit and the corresponding coordinates of each node are shown in the following formula. Figure 4 The unit with a length of l e is uniformly divided into n=2 j +m-2 parts, and the number of nodes of the unit is n+1. The actual coordinate values of the nodes in the unit are: x i ∈[x1,x n+1 ], 1≤i≤n+1.

[0071] The coordinate conversion formula is defined as: ξ=(x-x1) / l e , 0≤ξ≤1.

[0072] By substituting the actual coordinate values into the coordinate conversion formula, the mapping value ξ i of each node coordinate x i can be obtained.

[0073] ξ i =(x i -x1) / l e , 0≤ξ i ≤1, 1≤i≤n+1

[0074] When the m-order j-scale interval B-spline wavelet function is used as an interpolation function, the unknown field function u(ξ) can be expressed as:

[0075]

[0076] In the formula, a represents a coefficient vector.

[0077] Suppose the column vector of physical freedom is By substituting the mapping value ξ i of each node coordinate x i into the above formula, u e =R e a e .

[0078] In the formula, R e =[Φ T (ξ1)Φ T (ξ2)...Φ​T (ξ n+1 ] T , the unit conversion matrix is T e =(R e ) -1 .

[0079] Transforming u e =R e a e into a e =T e u e , we have

[0080]

[0081]

[0082] wherein represents a shape function matrix of a one-dimensional interval B-spline wavelet unit of m order and j scale.

[0083] 2. Discretize the finite element model of the plate structure to be measured by using the wavelet finite element method, divide the cross section Ω of the finite element model in the thickness direction into a finite number of units Ω e , and then solve the unit displacement of any point in the unit by using the interval B-spline wavelet unit shape function.

[0084] The harmonic displacement, stress and strain components of each point in each unit of the plate structure to be measured can be expressed as:

[0085] u x u y u z ] T ,σ x σ y σ z σ yz σ xz σ xy ] T ,ε x ε y ε z γ yz γ xz γ xy ] T

[0086] The constitutive relation of stress and strain is σ=Cε, wherein C is the elastic constant matrix of the plate structure to be measured. The relation of strain and displacement is:

[0087]

[0088] ​

[0089] Suppose the displacement field of ultrasonic guided wave along the propagation direction x is a simple harmonic wave, the displacement of any point in the unit is expressed by a spatial function, which is expressed as:

[0090]

[0091] When the scaling function of the m-order j-scale interval B-spline wavelet function is used as the interpolation function to construct the unit, there are n+1 nodes in each unit. The unit displacement of any point in the unit is expressed by the interval B-spline unit shape function. The unit displacement can be expressed as:

[0092]

[0093]

[0094] In the formula, Given by formula .

[0095] 3. According to the relationship between the unit displacement and the unit strain, the unit strain of each unit is obtained.

[0096] Substitute the unit displacement into the relationship between strain and displacement to obtain the unit strain of the unit:

[0097]

[0098]

[0099] 4. Substitute the unit displacement and the unit strain into the discrete form of the Hamilton formula to obtain the stiffness matrix and the mass matrix of each unit, respectively.

[0100] Let the number of units of the finite number of units Ω e obtained by discretizing the section Ω be n el , and introduce the discrete form of the Hamilton formula as:

[0101]

[0102] Substitute the unit displacement and the unit strain into the above formula to obtain:

[0103]

[0104]

[0105]

[0106]

[0107]

[0108] where t0 and t1 represent the start time and end time of the particle motion respectively, n el is the number of finite elements, e = 1, 2,..., n el , ε (e) represents the element strain, δ(ε (e)T ) represents the virtual strain, δ(u (e)T ) represents the virtual displacement, C e represents the material stiffness matrix, u (e) represents the element displacement, ρ e represents the material density of the plate structure to be measured, represents the second order inverse of the element displacement with respect to time t, (·) T represents the transpose of the matrix.

[0109] where, and are the three stiffness matrices corresponding to the element, m (e) is the mass matrix corresponding to the element, ρ e is the material density of the plate structure to be measured.

[0110] For the convenience of calculation, specific operator symbols Γ 0,0 , Γ 1,0 , Γ 1,0 and Γ 1,1 are introduced.

[0111]

[0112]

[0113] Using the specific operator, the stiffness matrix and the mass matrix are simplified. Take for example:

[0114]

[0115]

[0116]

[0117] The sub-matrix of Γ 1,1 can be calculated by the following formula:

[0118]

[0119] And m (e) can be calculated by the same method.

[0120] 5. The stiffness matrix and the mass matrix of each unit are integrated into a global stiffness matrix, which is substituted into the wave equation, and the wave equation is solved to obtain the dispersion curves of phase velocity and group velocity.

[0121] The stiffness matrix and the mass matrix are integrated into a global stiffness matrix and a global mass matrix by the following formula.

[0122]

[0123] In the formula, K1, K2 and K3 are the global stiffness matrices corresponding to and respectively, and M is the global mass matrix corresponding to m (e) .

[0124] The global stiffness matrix and the global mass matrix are substituted into the wave equation:

[0125] [K1+iκK2+κ 2 K3-ω 2 M] O U=0

[0126] In the formula, κ represents the wave number of the ultrasonic guided wave, i represents the imaginary symbol, K1, K2 and K3 are the global stiffness matrices corresponding to and respectively, M is the global mass matrix corresponding to m (e) , and U represents the displacement vector at the node, and subscript O is the total number of degrees of freedom.

[0127] Since the imaginary number i in the formula will bring great difficulty to the solution, a diagonal matrix T0 is introduced to eliminate the imaginary symbol i in the wave equation. If the element in the diagonal matrix T0 is related to u y , u z , then it is 1, and if it is related to u x , then it is i.

[0128]

[0129] In the formula, subscript O is the total number of degrees of freedom. Let i.e. Therefore, the wave equation can be expressed as:

[0130]

[0131] By solving the wave equation, 2M eigenvalues corresponding to the wave number k can be obtained for each frequency. Among all the eigenvalues, the number of real eigenvalues is the corresponding wave number. By using the relationship between the wave number and the frequency, and through certain numerical transformation, the dispersion curves of phase velocity and group velocity can be plotted.

[0132] The effect of the application is demonstrated by taking each isotropic homogeneous aluminum plate and each anisotropic composite material plate as examples.

[0133] The material parameters of each isotropic homogeneous aluminum plate are: thickness 1mm, elastic modulus 68.9Gpa, Poisson's ratio 0.33, and density 2690kg / m 3 . As shown in Figure 5 and Figure 6 , when the conventional finite element method is used for analysis, 10 units need to be constructed in the thickness direction, while after the application is used, only 3 units need to be constructed in the thickness direction, and the frequency dispersion curve obtained based on the application has a high degree of coincidence with the theoretical calculation result, but the calculation time is shortened from 66s to 25s, and on the basis of ensuring calculation accuracy, the calculation efficiency is greatly increased.

[0134] As shown in Figure 7 , the model of the anisotropic composite material plate, in order to further embody the advantage of the application that the calculation efficiency is high on the premise of ensuring calculation accuracy. When the conventional finite element method is used for analysis, 40 units are divided in the thickness direction, while after the application is used, 8 units are constructed in the thickness direction, and the total number of nodes is the same. As can be seen from Figure 8 and Figure 9 , compared with the frequency dispersion curve obtained by using the conventional finite element method, the group velocity and phase velocity frequency dispersion curves obtained after using the application have the same accuracy, and as can be seen from Figure 10 , after using the application, the calculation time is reduced from 476s to 226s, and the calculation efficiency is greatly increased.

Claims

1. A method for calculating the dispersion curve of an ultrasonic guided wave, characterized in that, Includes the following steps: S1. Establish a finite element model in the thickness direction of the plate structure to be tested, and divide the finite element model into elements; S2. Using the shape function constructed from the interval B-spline wavelet element, calculate the element displacement of each element in the finite element model; S3. The element displacement u calculated based on step S2 (e) The element strain ε of the corresponding element is calculated using the following formula. e ; S4. Calculate the stiffness matrix and mass matrix of the corresponding element based on the element displacement and element strain of each element. S5. Based on the stiffness matrix and mass matrix of each element, establish the wave equation and solve the wave equation to obtain the relationship between the wave number and frequency of the ultrasonic guided wave, and then draw the dispersion curve of the ultrasonic guided wave when it propagates in the structure under test.

2. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 1, characterized in that, The shape function constructed by the interval B-spline wavelet unit in step S2 is obtained by the following formula: In the formula, Let T be the shape function of the interval B-spline wavelet unit, Φ be the row vector of the scaling function of the interval B-spline wavelet function, and T be the shape function of the interval B-spline wavelet unit. e Let be the element transformation matrix of a one-dimensional interval B-spline wavelet element.

3. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 1, characterized in that, In step S2, the element displacement of the e-th element within the finite element method is calculated using the following formula: In the formula, u (e) For element displacement, Let q be the shape function of the interval B-spline wavelet unit. (e) e represents the nodal displacement. i(κx-ωt) For a simple harmonic wave, i is an imaginary number, κ is the wave number of the ultrasonic guided wave, x represents the propagation direction of the ultrasonic guided wave, ω is the frequency of the ultrasonic guided wave, and t is time.

4. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 1, characterized in that, In step S4, the stiffness matrix and mass matrix of each element are obtained by substituting the element displacement and element strain into the discrete form of the Hamiltonian formula; the discrete form of the Hamiltonian formula is as follows: In the formula, t0 and t1 represent the start and end times of the particle's motion, respectively, and n el Let e ​​be the number of elements in the finite element method, where e = 1, 2, ..., n. el , ε (e) δ(ε) represents the element strain. (e)T ) represents virtual strain, δ(u) (e)T ) represents virtual displacement, C e Represents the material stiffness matrix, u (e) ρ represents the element displacement. e This represents the material density of the structure of the plate under test. It represents the second reciprocal of the element displacement with respect to time t.

5. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 4, characterized in that, Substituting the element displacement and element strain into the discretized form of the Hamilton formula yields the following formula, which is then solved to obtain the stiffness matrix and mass matrix of the corresponding element. In the formula, δq (e)T Indicates virtual displacement. and For the three stiffness matrices within the corresponding element, m (e) Here is the mass matrix of the corresponding element, κ represents the wavenumber of the ultrasonic guided wave, ω represents the frequency of the ultrasonic guided wave, and q... (e) Represents the nodal displacement, where i is an imaginary number.

6. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 5, characterized in that, In step S5, when establishing the wave equation, a specific operator symbol Γ is introduced. 0,0 ,Γ 1,0 ,Γ 1,0 and Γ 1,1 Simplify the stiffness matrix and mass matrix; In the formula, T e Represents the element transformation matrix, l e Let Φ represent the length of the e-th cell, and let Φ represent the row vector of the scaling function of the interval B-spline wavelet function.

7. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 6, characterized in that, The stiffness matrix and mass matrix are integrated into a global stiffness matrix and a global mass matrix using the following formula, and then the wave equation is established using the global stiffness matrix and the global mass matrix. In the formula, K1, K2, and K3 are respectively the corresponding and The overall stiffness matrix, M is the corresponding m (e) The overall quality matrix.

8. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 7, characterized in that, The wave equation is: [K1+iκK2+κ 2 K3-oh 2 M] O U=0 In the formula, κ represents the wavenumber of the ultrasonic guided wave, i represents the imaginary number, and K1, K2, and K3 are the corresponding... and The overall stiffness matrix, M is the corresponding m (e) The overall mass matrix is ​​given by , where U represents the displacement vector at the node and the subscript O is the total number of degrees of freedom.

9. The method for calculating the dispersion curve of ultrasonic guided waves according to claim 8, characterized in that, When solving the wave equation, a diagonal matrix T0 is introduced to eliminate the imaginary symbol i in the wave equation. The diagonal matrix T0 is: In the formula, the subscript O represents the total number of degrees of freedom.

Images