A meshless buckling analysis method for variable stiffness plates with delamination based on delamination theory

Abstract

The present invention discloses a meshless buckling analysis method for a delaminated variable stiffness plate based on layered theory, comprising: obtaining geometric and material parameter information of a delaminated variable stiffness fiber composite laminate; generating node and integration point information based on the parameter information combined with numerical layer number information at the delaminated location, calculating and storing meshless shape functions and their derivative information of the nodes and integration points, and completing meshless spatial discretization of the problem domain; constructing an overall stiffness matrix and a geometric stiffness matrix based on layered approximation theory, calculating the buckling eigenvalues ​​and buckling eigenvectors of the delaminated variable stiffness fiber composite laminate using a direct method to impose essential boundary conditions; analyzing the penetration behavior of the delaminated variable stiffness fiber composite laminate based on the buckling eigenvalues ​​and buckling eigenvectors. The method of the present invention is applicable to the buckling analysis of curved fiber composite laminates with arbitrary delamination; it can eliminate the mutual intrusion of the upper and lower parts of the delaminated area, and the results are more reasonable and accurate.

Description

Technical Field

[0001] The present invention relates to the technical field of variable stiffness fiber composite materials, and in particular to a gridless buckling analysis method for variable stiffness plates containing delamination based on delamination theory. Background Art

[0002] Resin-based fiber composites have the advantages of high specific strength, high specific stiffness, corrosion resistance, and easy processing and molding. They are widely used in a variety of fields such as aerospace, automobiles, ships, and wind power. Traditional constant-stiffness fiber composite laminates use straight fibers, while variable-stiffness fiber composites use curved fiber layups, which provide more design space and can achieve better mechanical properties. Moreover, with the development of advanced automatic fiber placement technology in recent years, curved fiber placement has moved from theoretical ideas to practical preparation. It is now possible to change the placement angle of composite materials by placing fiber bundles in a curved manner, achieving a more reasonable force transmission path.

[0003] Resin-based fiber composite laminates are typically used in engineering applications as plate and shell structures. As the primary failure mode of composite plate and shell structures, the calculation of buckling strength is particularly important in design. Calculation of variable-stiffness laminates is even more complex due to the curved fiber paths. Furthermore, due to the process and structural characteristics of fiber composite laminates, the interlaminar strength, which primarily depends on the resin substrate and the bond between the substrate and the fibers, is far lower than the intralaminar fiber strength. Delamination failure is one of the primary failure modes of fiber composite laminates. Currently, no effective method has been developed for calculating the buckling bearing capacity of variable-stiffness fiber composite laminates with delamination. Summary of the Invention

[0004] The purpose of this section is to summarize some aspects of the embodiments of the present invention and briefly introduce some preferred embodiments. Some simplifications or omissions may be made in this section and the abstract and title of this application to avoid obscuring the purpose of this section, the abstract and the title of the invention, and such simplifications or omissions should not be used to limit the scope of the present invention.

[0005] In view of the above existing problems, the present invention is proposed.

[0006] Therefore, the present invention provides a meshless buckling analysis method for variable stiffness plates with delamination based on delamination theory to solve the existing problem that it is impossible to correctly analyze and predict the buckling bearing capacity and buckling mode of variable stiffness fiber composite laminates with delamination.

[0007] In order to solve the above technical problems, the present invention provides the following technical solutions:

[0008] The embodiment of the present invention provides a meshless buckling analysis method for a delamination-containing variable stiffness plate based on delamination theory, including:

[0009] Obtain the geometric and material parameter information of variable stiffness fiber composite laminates with delamination;

[0010] Based on the parameter information and the numerical layer number information at the delamination point, node and integration point information is generated, and meshless shape functions and derivative information of the nodes and integration points are calculated and stored to complete meshless spatial discretization of the problem domain;

[0011] Based on the layered approximation theory, the global stiffness matrix and the geometric stiffness matrix are constructed. The buckling eigenvalues ​​and eigenvectors of variable stiffness fiber composite laminates with delamination are calculated by applying essential boundary conditions through the direct method.

[0012] The penetration behavior of variable stiffness fiber composite laminates with delamination is analyzed based on the buckling eigenvalue and buckling eigenvector.

[0013] As a preferred solution of the meshless buckling analysis method of delamination-containing variable stiffness plates based on the delamination theory of the present invention, it also includes:

[0014] If there is no element with a relative displacement smaller than the first threshold at the delamination interface, it is determined that there is no penetration behavior;

[0015] If there is an element whose relative displacement at the delamination interface is smaller than the first threshold, it is determined that penetration behavior exists;

[0016] When there is penetration behavior and the penetration displacement is not greater than the allowable penetration threshold, the penetration situation is ignored;

[0017] When penetration behavior exists, if the penetration displacement is greater than the allowable penetration threshold, the iteration parameters are determined, the spring additional stiffness matrix is ​​calculated, and the overall stiffness matrix is ​​updated until there is no penetration behavior.

[0018] As a preferred solution of the meshless buckling analysis method of delaminated variable stiffness plates based on the delamination theory described in the present invention, the calculation of the spring additional stiffness matrix includes: the stiffness of the virtual spring at any node i of the delaminated area is expressed as:

[0019]

[0020]

[0021] in, is the delamination interface point x in the numerical solution of this iteration i The relative displacement of the interface at The diagonal line of the overall stiffness matrix corresponds to Item, Point x after stiffness change i The relative displacement of the interface at , r is the iteration parameter and r=0.00001×3 m-1, m is the number of iterations, and i is any node with penetration behavior.

[0022] As a preferred solution of the meshless buckling analysis method of delamination-containing variable stiffness plates based on the delamination theory described in the present invention, the overall stiffness matrix is ​​expressed as:

[0023]

[0024] Among them, U, are the displacement parameters of the laminate delamination region and the delamination region, respectively, K b , K s , K t are the non-delaminated bending matrix, shear matrix, and coupled stiffness matrix, respectively. are the bending matrix, shear matrix, and coupling stiffness matrix of delamination and non-delamination coupling, are the bending matrix, shear matrix, and coupling stiffness matrix of the delamination, respectively, and K is the overall stiffness matrix.

[0025] As a preferred solution of the meshless buckling analysis method of delamination-containing variable stiffness plates based on the layered theory of the present invention, the geometric stiffness matrix is ​​calculated based on the work done by the pre-buckling loads of NP numerical layers, including:

[0026] The geometric stiffness matrix K G Expressed as:

[0027]

[0028] in, is the deflection of the nth numerical layer, is the in-plane horizontal axial compressive prebuckling load of the nth numerical layer, is the axial compressive prebuckling load in the direction perpendicular to the plane coordinate of the nth numerical layer, is the in-plane shear prebuckling load, W1 is the geometric stiffness matrix of the nth numerical layer, W4 is the geometric stiffness matrix of the n+1th numerical layer, W2 is the first coupling stiffness matrix, and W3 is the second coupling stiffness matrix.

[0029] As a preferred solution of the meshless buckling analysis method for variable stiffness plates with delamination based on the delamination theory of the present invention, assuming that delamination exists at both the nth and n+1th numerical interfaces, the derivative of the numerical interface deflection can be expressed by the meshless shape function as follows:

[0030]

[0031]

[0032] Where ne is the number of functions, is the displacement parameter corresponding to the meshless node, W n,x is the derivative of the deflection of the nth numerical layer in the x direction, W n,y is the derivative of the deflection of the nth numerical layer in the y direction, W n+1,x is the derivative of the deflection of the n+1th numerical layer in the x direction, W n+1,y is the derivative of the deflection of the n+1th numerical layer in the y direction, is the shape function matrix of the i-th node, N i,x is the derivative of the shape function of the node at position i with respect to x, N i,y is the derivative of the shape function of the node at position i with respect to y, is the displacement parameter in the u direction of the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, is the displacement parameter in the v direction of the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, is the displacement parameter in the w direction at the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, and the superscript “^” indicates delamination.

[0033] As a preferred solution of the meshless buckling analysis method of delamination-containing variable stiffness plates based on the delamination theory of the present invention, it also includes:

[0034] If there is no delamination at the nth or n+1th numerical interface, then the corresponding or The item does not exist;

[0035] The direct method is used to impose the essential boundary conditions. At the boundary node i, let The calculated buckling eigenvalue and buckling mode are expressed as:

[0036] (K-λK G )U=0

[0037] Among them, λ is the eigenvalue and U is the eigenvector.

[0038] As a preferred solution of the meshless buckling analysis method of the delamination-containing variable stiffness plate based on the delamination theory of the present invention, the geometric and material parameter information of the delamination-containing variable stiffness fiber composite material laminate is obtained, which at least includes: the plate length L x , board width L y , plate thickness H, number of layers NL, layer thickness t i , peeling quantity, peeling position h / H, peeling width L dx / L x , L dy / Ly, elastic modulus, shear modulus, Poisson's ratio.

[0039] As a preferred solution of the meshless buckling analysis method of delamination-containing variable stiffness plates based on the delamination theory of the present invention, it also includes:

[0040] Input fiber routing and layup design

[0041] According to the board length L x and board width L y , as well as the number of numerical layers, layout nodes and integration points.

[0042] As a preferred solution of the meshless buckling analysis method for delamination-containing variable stiffness plates based on the delamination theory described in the present invention, the method includes: calculating and storing meshless shape functions and their derivative information of nodes and integration points to complete meshless spatial discretization of the problem domain, including: determining the problem domain according to the geometric dimensions of the plate, and regularly arranging an appropriate number of nodes within the problem domain;

[0043] The radial base point interpolation method is used to calculate the shape function of each node and approximate the displacement field u h (x) is expressed as:

[0044] u h (x) = {rp}G -1 U=ΦU

[0045] Φ={rp}G -1

[0046] Among them, Φ is the shape function, r and p are radial basis functions and polynomial basis functions, and G is the basis function coefficient.

[0047] Compared with existing technologies, the present invention offers the following advantages: The method is applicable to the buckling analysis of curved fiber composite laminates with arbitrary delaminations; it eliminates the mutual intrusion between the upper and lower parts of the delamination region, resulting in more reasonable and accurate results; and it more accurately characterizes the variable stiffness characteristics of the structure, handling arbitrary plate and shell geometries and arbitrary delamination geometries without requiring quadratic approximation of the field function, resulting in improved computational efficiency. Furthermore, the method can analyze global, mixed, and local buckling, while also eliminating the potential penetration between the upper and lower layers of the delamination region during calculations. BRIEF DESCRIPTION OF THE DRAWINGS

[0048] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the following briefly introduces the drawings required for describing the embodiments. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be derived from these drawings without inventive effort. Among them:

[0049] Figure 1 A schematic flow chart of a meshless buckling analysis method for a delamination-containing variable stiffness plate based on delamination theory according to an embodiment of the present invention;

[0050] Figure 2 Schematic diagram of Model A of a meshless buckling analysis method for delamination-containing variable stiffness plates based on delamination theory according to one embodiment of the present invention;

[0051] Figure 3 Schematic diagram of Model B of the meshless buckling analysis method for delamination-containing variable stiffness plates based on delamination theory according to one embodiment of the present invention;

[0052] Figure 4 A linearly varying fiber path diagram of a meshless buckling analysis method for a delamination-containing variable stiffness plate based on delamination theory according to an embodiment of the present invention;

[0053] Figure 5 Schematic diagram of the loading mode of Model A of the meshless buckling analysis method for delamination-containing variable stiffness plates based on delamination theory according to one embodiment of the present invention;

[0054] Figure 6 Schematic diagram of the loading mode of Model B of the meshless buckling analysis method for delamination-containing variable stiffness plates based on delamination theory according to one embodiment of the present invention;

[0055] Figure 7 Buckling mode diagram of a delamination vertical position and width of 0.1 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0056] Figure 8 Buckling modal diagram of a delamination vertical position of 0.1 and a width of 0.5 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0057] Figure 9 Buckling modal diagram of a delamination vertical position of 0.1 and a width of 0.9 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0058] Figure 10 Buckling modal diagram of a delamination vertical position of 0.3 and a width of 0.1 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0059] Figure 11 Buckling modal diagram of a delamination vertical position of 0.3 and a width of 0.5 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0060] Figure 12Buckling modal diagram of a delamination vertical position of 0.3 and a width of 0.9 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0061] Figure 13 Buckling modal diagram of a delamination vertical position of 0.5 and a width of 0.1 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0062] Figure 14 Buckling modal diagram of a delamination vertical position of 0.5 and a width of 0.5 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0063] Figure 15 Buckling modal diagram of a delamination vertical position of 0.5 and a width of 0.9 according to a meshless buckling analysis method for a variable stiffness plate with delamination based on delamination theory according to an embodiment of the present invention;

[0064] Figure 16 Buckling modal diagram of Model B of the meshless buckling analysis method for delaminated variable stiffness plates based on delamination theory according to one embodiment of the present invention, using an additional spring to eliminate penetration phenomena;

[0065] Figure 17 This is a buckling modal diagram of Model B of the meshless buckling analysis method for delaminated variable stiffness plates based on delamination theory according to one embodiment of the present invention, in which no additional spring is used to eliminate penetration phenomena. DETAILED DESCRIPTION

[0066] To make the above-mentioned objects, features, and advantages of the present invention more clearly understood, the following detailed description of the specific embodiments of the present invention is given in conjunction with the accompanying drawings. It is obvious that the described embodiments are only part of the embodiments of the present invention, not all of them. Based on the embodiments of the present invention, all other embodiments obtained by ordinary persons in this field without creative work should fall within the scope of protection of the present invention.

[0067] In the following description, many specific details are set forth to facilitate a full understanding of the present invention. However, the present invention may also be implemented in other ways different from those described herein. Those skilled in the art may make similar generalizations without violating the connotation of the present invention. Therefore, the present invention is not limited to the specific embodiments disclosed below.

[0068] Secondly, the term "one embodiment" or "embodiment" herein refers to a specific feature, structure, or characteristic that may be included in at least one implementation of the present invention. The phrase "in one embodiment" appearing in various places throughout this specification does not necessarily refer to the same embodiment, nor does it refer to a separate or selective embodiment that is mutually exclusive of other embodiments.

[0069] The present invention is described in detail with reference to schematic diagrams. For ease of illustration, cross-sectional views of device structures may be partially enlarged and not to scale when describing embodiments of the present invention. Furthermore, the schematic diagrams are merely illustrative and should not limit the scope of the present invention. Furthermore, in actual production, the three-dimensional dimensions of length, width, and depth should be included.

[0070] In the description of the present invention, it should be noted that the terms "upper, lower, inner, and outer" and other references to orientations or positional relationships are based on the orientations or positional relationships shown in the accompanying drawings and are intended solely to facilitate and simplify the description of the present invention. They are not intended to indicate or imply that the devices or components referred to must have, be constructed, or operate in a specific orientation, and therefore should not be construed as limitations on the present invention. Furthermore, the terms "first, second, or third" are used for descriptive purposes only and should not be construed as indicating or implying relative importance.

[0071] In this disclosure, unless otherwise specified or limited, the terms "mounted," "connected," and "connected" should be interpreted broadly. For example, they may refer to fixed, removable, or integral connections. They may also refer to mechanical, electrical, or direct connections, indirect connections through an intermediary, or internal communication between two components. Those skilled in the art will understand the specific meanings of these terms in this disclosure.

[0072] Example 1

[0073] Reference Figure 1 , is an embodiment of the present invention, which provides a meshless buckling analysis method for variable stiffness plates with delamination based on delamination theory, comprising:

[0074] S1: Obtaining geometric and material parameter information of a delamination-variable stiffness fiber composite laminate;

[0075] Furthermore, the geometric and material parameter information of the delamination variable stiffness fiber composite laminate includes at least: the plate length L x , board width L y , plate thickness H, number of layers NL, layer thickness t i , peeling quantity, peeling position h / H, peeling width L dx / L x , L dy / Ly, elastic modulus, shear modulus, Poisson's ratio;

[0076] Furthermore, the fiber path and ply design are input. In the embodiment of the present invention, the fiber path is predetermined by the function θ(x, y) to analyze any fiber path and ply design. The linear fiber path is expressed as:

[0077]

[0078] Where x and y are the coordinates of the calculation point, φ is the rotation angle of the coordinate axis, d0 represents the linear change period of the fiber, T0 and T1 are the angles between the fiber and the x′ axis at the beginning and end of the period, respectively, and θ is the fiber deflection angle.

[0079] Specifically, the fiber path is recorded as [φ <t0t1>].

[0080] S2: Based on the parameter information and the numerical layer number information at the delamination point, the node and integration point information is generated, the meshless shape functions and their derivative information of the nodes and integration points are calculated and stored, and the meshless spatial discretization of the problem domain is completed;

[0081] Furthermore, the embodiment of the present invention adopts Gaussian integral according to the plate length L x and board width L y , as well as the number of numerical layers, arrangement nodes and integration points, the number of numerical layers used in the layered approximation theory can be greater than, equal to, or less than the number of plies.

[0082] Furthermore, the meshless shape functions and their derivatives of nodes and Gaussian points are calculated and stored for standby adjustment;

[0083] Specifically, the problem domain is determined according to the geometric dimensions of the plate, and appropriate nodes are arranged regularly within the problem domain. The shape function Φ of each node is calculated using RPIM, and the approximate function is set as

[0084]

[0085] Where n is the number of nodes in the support domain of the computational point, m is the number of terms in the polynomial basis function, and r i (x) and p j (x) are radial basis functions and polynomial basis functions, respectively, a i (x) and b j (x) is the corresponding coefficient.

[0086] Let the nodes in the computation point support domain satisfy the above approximate formula, and we can get:

[0087]

[0088] Where u is the column vector of the field function, R and P are matrices consisting of the row vectors of the radial basis function and the row vectors of the polynomial basis function, respectively; thus, we can obtain:

[0089] u h (x) = {rp}G -1 U=ΦU

[0090] Φ={rp}G -1

[0091] S3: Based on the layered approximation theory, the global stiffness matrix and the geometric stiffness matrix are constructed. The buckling eigenvalues ​​and buckling eigenvectors of the variable stiffness fiber composite laminate with delamination are calculated by applying the essential boundary conditions through the direct method.

[0092] Furthermore, the global stiffness matrix is ​​calculated based on the matrix form of the strain energy variation of the lower laminate, which is expressed as:

[0093]

[0094] Among them, U, are the displacement parameters of the laminate delamination region and the delamination region, respectively, K b , K s , K t are the non-delaminated bending matrix, shear matrix, and coupled stiffness matrix, respectively. are the bending matrix, shear matrix, and coupling stiffness matrix of delamination and non-delamination coupling, are the bending matrix, shear matrix, and coupling stiffness matrix of the delamination, respectively, and K is the overall stiffness matrix.

[0095] Specifically, the non-delaminated bending matrix is ​​expressed as:

[0096]

[0097]

[0098] Where T represents the matrix transpose, NP is the number of numerical layers, M is the numerical layer plane and M=NP+1, are the ordinates of the lower and upper surfaces of the kth numerical layer, is a non-delamination elastic matrix, is the non-delamination bending matrix at the I, J numerical plane, H I , H J is the interpolation function at the I, J numerical plane, is the bending elastic coefficient of the k-value layer, B b is the strain matrix of the non-delamination bending matrix, D b is the strain matrix, H I (z) is the global interpolation function of the displacement component along the thickness direction, Ψ I (z) is a step function representing delamination.

[0099] The bending matrix of delamination and non-delamination coupling is expressed as:

[0100]

[0101]

[0102] in, is the strain matrix of the delaminated and non-delaminated coupled bending matrix, is the bending stiffness matrix of the I and J numerical plane layers coupled with delamination and non-delamination.

[0103] The bending matrix of the delamination is expressed as:

[0104]

[0105] in, is the bending stiffness moment of the delaminated layer I and J.

[0106] The non-delamination shear matrix is ​​expressed as:

[0107]

[0108]

[0109] Among them, B s is the strain matrix of the non-delamination bending-shear coupling stiffness moment, is the bending-shear coupling stiffness moment of the numerical plane layer of the non-delaminated layer I and J.

[0110] The shear stiffness matrix of delamination and non-delamination coupling is expressed as:

[0111]

[0112]

[0113] in, The strain matrix of the bending-shear coupling stiffness moment of delamination and non-delamination coupling, is the bending-shear coupling stiffness moment of the I and J numerical plane layers of the non-delamination and delamination coupling, are the strain matrices of the i-th shape function respectively.

[0114] The shear stiffness matrix of delamination is expressed as:

[0115]

[0116] in, is the bending-shear coupling stiffness moment of the delaminated layers I and J.

[0117] The non-delaminated shear stiffness matrix is ​​expressed as:

[0118]

[0119] Among them, B t is the strain matrix of the non-delamination shear stiffness matrix,

[0120] The shear stiffness matrix of delamination and non-delamination coupling is expressed as:

[0121]

[0122]

[0123] in, is the strain matrix of the shear stiffness matrix of the delamination and non-delamination coupling,

[0124] The shear stiffness matrix of delamination is expressed as:

[0125]

[0126]

[0127] in, is the strain matrix of the i-th shape function.

[0128] Furthermore, if the global interpolation function H I (z) is the Lagrange interpolation function, which can be expressed as follows

[0129]

[0130]

[0131]

[0132]

[0133] Among them, h k is the thickness of the kth numerical layer, is the thickness coordinate of the lower surface of the kth numerical layer, and z is the vertical coordinate.

[0134] Furthermore, the step function of delamination is expressed as:

[0135]

[0136] Among them, z I The vertical coordinate represents the first delamination.

[0137] is the strain matrix of the i-th shape function, expressed as follows:

[0138]

[0139]

[0140]

[0141] Among them, φ i 、φ i,x and φ i,y are the meshless shape function and its derivatives in the x and y directions, respectively.

[0142] and They are the elastic coefficients of the in-plane bending, transverse shear, in-plane bending and transverse expansion coupling parts of the k-th numerical layer in the global coordinate system, which are expressed as follows

[0143]

[0144] The stiffness coefficient matrix of the kth layer of the laminate in the global coordinate system can be obtained by the stiffness matrix of the material coordinate system through the following coordinate transformation, which is expressed as:

[0145]

[0146]

[0147] Here, θ is the fiber angle.

[0148] Furthermore, based on the work done by the prebuckling load of the NP numerical layers, the geometric stiffness matrix is ​​calculated including:

[0149] Geometric stiffness matrix K G Expressed as:

[0150]

[0151] in, is the deflection of the nth numerical layer, is the in-plane horizontal axial compressive prebuckling load of the nth numerical layer, is the axial compressive prebuckling load in the direction perpendicular to the plane coordinate of the nth numerical layer, is the in-plane shear prebuckling load, W1 is the geometric stiffness matrix of the nth numerical layer, W4 is the geometric stiffness matrix of the n+1th numerical layer, W2 is the first coupling stiffness matrix, and W3 is the second coupling stiffness matrix.

[0152] Specifically, the geometric stiffness matrix of the nth numerical layer is expressed as:

[0153]

[0154] Among them, W n,x ,W n,y is the derivative of the deflection of the nth numerical layer in the x and y directions.

[0155] The first coupling stiffness matrix and the second coupling stiffness matrix are expressed as:

[0156]

[0157] The geometric stiffness matrix of the n+1th numerical layer is expressed as:

[0158]

[0159] Furthermore, assuming that delamination exists at both the nth and n+1th numerical interfaces, the derivative of the numerical interface deflection can be expressed by the meshless shape function as follows:

[0160]

[0161]

[0162] Where ne is the number of shape functions, is the corresponding displacement parameter, W n,x is the derivative of the deflection of the nth numerical layer in the x direction, W n,y is the derivative of the deflection of the nth numerical layer in the y direction, W n+1,x is the derivative of the deflection of the n+1th numerical layer in the x direction, W n+1,y is the derivative of the deflection of the n+1th numerical layer in the y direction, is the shape function matrix of the i-th node, N i,x is the derivative of the shape function of the node at position i with respect to x, N i,y is the derivative of the shape function of the node at position i with respect to y, is the displacement parameter in the u direction of the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, is the displacement parameter in the v direction of the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, is the displacement parameter in the w direction at the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, and the superscript “^” indicates delamination.

[0163] Furthermore, if there is no delamination at the nth or n+1th numerical interface, then the corresponding or The item does not exist.

[0164] The direct method is used to impose the essential boundary conditions. At the boundary node i, let For the system equation (K-λK G )U=0 and solve to obtain the buckling eigenvalue and buckling mode.

[0165] S4: Analysis of the penetration behavior of variable stiffness fiber composite laminates with delamination based on buckling eigenvalues ​​and buckling eigenvectors.

[0166] Furthermore, it also includes:

[0167] If there is no element with a relative displacement smaller than the first threshold at the delamination interface, it is determined that there is no penetration behavior;

[0168] If there is an element whose relative displacement at the delamination interface is smaller than the first threshold, it is determined that penetration behavior exists;

[0169] When there is penetration behavior and the penetration displacement is not greater than the allowable penetration threshold, the penetration situation is ignored;

[0170] When penetration behavior exists, if the penetration displacement is greater than the allowable penetration threshold, the iteration parameters are determined, the spring additional stiffness matrix is ​​calculated, and the overall stiffness matrix is ​​updated until there is no penetration behavior.

[0171] It should be noted that, in the embodiment of the present invention, the value of the first threshold is 0, and the permissible penetration threshold range is [10 -4 ,10 -6 ], the specific value of the allowable penetration threshold can be selected according to the model size and the specific problem.

[0172] Specifically, the calculation of the spring additional stiffness matrix includes: the stiffness of the virtual spring at any node of the delamination area is expressed as:

[0173]

[0174]

[0175] in, is the delamination interface point x in the numerical solution of this iteration i The relative displacement of the interface at The diagonal line of the overall stiffness matrix corresponds to Item, Point x after stiffness change i The relative displacement of the interface at , r is the iteration parameter and r=0.00001×3 m-1 , m is the number of iterations, and i is any node with penetration behavior.

[0176] It should be noted that setting r to a parameter that increases exponentially with the number of iterations is intended to minimize the effect of the additional stiffness on the displacements of nodes at other locations when the number of iterations is low, while facilitating faster convergence when the number of iterations is high. The choice of parameter r determines the speed of iterative convergence and the accuracy of the iterative results.

[0177] Example 2

[0178] Reference Figures 2 to 17 , which is an embodiment of the present invention. Different from the first embodiment, this embodiment provides a specific example of the gridless buckling analysis method of delamination-containing variable stiffness plates based on the delamination theory, verifying the beneficial effects of our invention.

[0179] This embodiment provides two different delamination situations, namely model A and model B.

[0180] like Figure 2 The figure shows the geometric dimensions of model A. Model A has a longitudinal through-delamination, where L x =4m, L y =4m, H = 0.01m. The boundary conditions are simple support on the opposite side and free on the opposite side; the thickness of each layer is t i =0.001m. Delamination position h / H=0.1, 0.3, 0.5; Delamination width L d / L x =0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9;

[0181] like Figure 3 The figure shows the geometric dimensions of model B. There is a rectangular delamination area in the middle of model B, where: L x =L y =0.81408m, H=1.016mm. Its delamination position h / H=0.25, delamination width L dx / L x =L dy / L y =0.6. The boundary conditions are simply supported on the opposite side and free on the opposite side;

[0182] In the material coordinate system, the elastic modulus of the laminate is: E 11 =181Gpa, E 22 =E 33 =10.273GPa, shear modulus G 12 =G 23 =G 13 =7.170GPa, Poisson's ratio ν 12 =0.28. Input fiber path definition function θ(x,y): Figure 4 As shown, the fiber laying state of model A laminate is [±0<0|0>] 10 ,like Figure 5 As shown, the fiber orientation of model B is [±90<0|75>] 2s .

[0183] According to the board length L x , width L y , Model A is divided into 40 equal-sized integration units, and Model B is divided into 64 integration units. Both are calculated using two numerical layers. 4-point Gaussian integration is used to calculate the integration point coordinates, weights, and Jacobian ratios; the RPIM calculation shape function parameters q = 1.5, α c =1, dimensionless support domain size d max =5, regularly arrange 101×11 nodes in the problem domain, calculate the shape functions of the nodes and store them.

[0184] The stiffness matrix K of models A and B is calculated based on the matrix form of the strain energy variation of the laminate, and the work done by the longitudinal prebuckling load is integrated for the two numerical layers to obtain the geometric stiffness matrix K of models A and B. G ; Use the direct method to impose the essential boundary conditions, at the boundary node i, let For the system equation (K-λK G )U=0 and solve to obtain the buckling eigenvalue and buckling mode.

[0185] Determine whether there is penetration behavior, that is, the relative displacement at the delamination interface Is there an element with a displacement smaller than the preset allowable penetration displacement - ε? If not, it is determined that there is no penetration behavior and the process ends. The allowable penetration displacement ε is 0.00001. If there is penetration behavior, calculate the spring stiffness K. * , update the global stiffness matrix K, solve the buckling eigenvalues ​​and buckling modes, and end when there is no penetration. Based on the buckling eigenvectors, draw the buckling mode diagram.

[0186] In this example, the nodes are arranged in a regular 101×11 pattern. FEM results were obtained using Abaqus software using continuous shell element modeling, employing 100×100 elements. Under Model A, a comparison of normalized buckling loads was performed using the radial basis meshless layered method (RPIM-LW) of the present invention, a one-dimensional model established by GJ Simitses et al., and a finite element method (FEM) model. The buckling performance of the laminate was characterized using the defined nondimensional buckling load, expressed as:

[0187]

[0188] in, is the buckling load, L x is the board length, E 11 is the elastic modulus in the x direction, and H is the plate thickness.

[0189] The comparison results are shown in Table 1:

[0190] Table 1 Normalized buckling load of model A laminate with delamination

[0191]

[0192]

[0193] Combine Figures 7 to 15 As can be seen from Table 1, the buckling loads calculated by the three methods for each model are basically consistent, proving the effectiveness of the method of the present invention. Figure 7 As shown, it indicates that the delamination is at the vertical position of h / H=0.1, and the width of the delamination is L d / L x =0.1, i.e. longitudinal full-through delamination with a delamination width of 0.4 m.

[0194] Based on the method of the present invention under the conditions of model B, the fiber orientation is calculated to be [±90<0|75>] 2s The normalized buckling loads of the VSC laminates with and without the additional spring to remove penetration are shown in Table 2:

[0195] Table 2 Normalized buckling load of model B laminate with delamination

[0196] Are there any constraints? Fiber orientation RPIM-LW Unrestrained buckling <![CDATA[[±90<0|75>] 2s ]]> 0.8038 Considering contact buckling <![CDATA[[±90<0|75>] 2s ]]> 1.1454

[0197] Combine Figures 16-17 From Table 2, it can be seen that the present invention is applicable to the buckling analysis of curved fiber composite laminates with arbitrary delamination; the mutual intrusion of the upper and lower parts of the delamination area can be eliminated, and the results are more reasonable and accurate.

[0198] It should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit the present invention. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present invention may be modified or replaced by equivalents without departing from the spirit and scope of the technical solutions of the present invention, which should all be included in the scope of the claims of the present invention.

Claims

1. A meshless buckling analysis method for variable stiffness plates with delamination based on delamination theory, characterized by: include: Obtain the geometric and material parameter information of variable stiffness fiber composite laminates with delamination; Based on the parameter information and the numerical layer number information at the delamination point, node and integration point information is generated, and meshless shape functions and derivative information of the nodes and integration points are calculated and stored to complete meshless spatial discretization of the problem domain; Based on the layered approximation theory, the global stiffness matrix and geometric stiffness matrix are constructed. The buckling eigenvalues ​​and buckling eigenvectors of variable stiffness fiber composite laminates with delamination are calculated by applying essential boundary conditions through the direct method. The global stiffness matrix is ​​expressed as: Among them, U, are the displacement parameters of the laminate delamination region and the delamination region, respectively, K b , K s , K t are the non-delaminated bending matrix, shear matrix, and coupled stiffness matrix, respectively. are the bending matrix, shear matrix, and coupling stiffness matrix of delamination and non-delamination coupling, are the bending matrix, shear matrix, and coupled stiffness matrix of the delamination, respectively, and K is the overall stiffness matrix; According to the work done by the prebuckling load of NP numerical layers, the geometric stiffness matrix is ​​calculated including: The geometric stiffness matrix K G Expressed as: in, is the deflection of the nth numerical layer, is the in-plane horizontal axial compressive prebuckling load of the nth numerical layer, is the axial compressive prebuckling load in the direction perpendicular to the plane coordinate of the nth numerical layer, is the in-plane shear prebuckling load, W1 is the geometric stiffness matrix of the nth numerical layer, W4 is the geometric stiffness matrix of the n+1th numerical layer, W2 is the first coupling stiffness matrix, and W3 is the second coupling stiffness matrix; Assuming that delamination exists at both the nth and n+1th numerical interfaces, the derivative of the numerical interface deflection can be expressed by the meshless shape function as follows: Where ne is the number of functions, is the displacement parameter corresponding to the meshless node, W n,x is the derivative of the deflection of the nth numerical layer in the x direction, W n,y is the derivative of the deflection of the nth numerical layer in the y direction, W n+1,x is the derivative of the deflection of the n+1th numerical layer in the x direction, W n+1,y is the derivative of the deflection of the n+1th numerical layer in the y direction, is the shape function matrix of the i-th node, N i,x is the derivative of the shape function of the node at position i with respect to x, N i,y is the derivative of the shape function of the node at position i with respect to y, is the displacement parameter in the u direction of the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, is the displacement parameter in the v direction of the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, is the displacement parameter in the w direction at the non-delamination point corresponding to the i-th shape function of the n+1-th numerical layer, and the superscript "^" indicates delamination; If there is no delamination at the nth or n+1th numerical interface, then the corresponding or The item does not exist; The direct method is used to impose the essential boundary conditions. At the boundary node i, let The calculated buckling eigenvalue and buckling mode are expressed as: (K-λK G )U=0 Among them, λ is the eigenvalue, U is the eigenvector; The penetration behavior of variable stiffness fiber composite laminates with delamination is analyzed based on the buckling eigenvalue and buckling eigenvector.

2. The meshless buckling analysis method for variable stiffness plates with delamination based on delamination theory according to claim 1 is characterized in that: Also includes: If there is no element with a relative displacement smaller than the first threshold at the delamination interface, it is determined that there is no penetration behavior; If there is an element whose relative displacement at the delamination interface is smaller than the first threshold, it is determined that penetration behavior exists; When there is penetration behavior and the penetration displacement is not greater than the allowable penetration threshold, the penetration situation is ignored; When penetration behavior exists, if the penetration displacement is greater than the allowable penetration threshold, the iteration parameters are determined, the spring additional stiffness matrix is ​​calculated, and the overall stiffness matrix is ​​updated until there is no penetration behavior.

3. The meshless buckling analysis method for variable stiffness plates with delamination based on delamination theory according to claim 1 or 2, characterized in that: The calculation of the spring additional stiffness matrix includes: the stiffness of the virtual spring at any node i in the delamination area is expressed as: in, is the delamination interface point x in the numerical solution of this iteration i The relative displacement of the interface at The diagonal line of the overall stiffness matrix corresponds to Item, Point x after stiffness change i The relative displacement of the interface at , r is the iteration parameter and r=0.00001×3 m-1 , m is the number of iterations, and i is any node with penetration behavior.

4. The meshless buckling analysis method for variable stiffness plates with delamination based on delamination theory according to claim 3 is characterized in that: The geometric and material parameter information of the delamination variable stiffness fiber composite laminate includes at least: plate length L x , board width L y , plate thickness H, number of layers NL, layer thickness t i , peeling quantity, peeling position h / H, peeling width L dx / L x , L dy / Ly , elastic modulus, shear modulus, Poisson's ratio.

5. The meshless buckling analysis method for variable stiffness plates with delamination based on layered theory according to claim 4 is characterized in that: Also includes: Input fiber routing and layup design According to the board length L x and board width L y , as well as the number of numerical layers, layout nodes and integration points.

6. The meshless buckling analysis method for variable stiffness plates with delamination based on layered theory according to claim 5, characterized in that: Calculate and store the meshless shape functions and their derivatives of nodes and integration points to complete the meshless spatial discretization of the problem domain, including: determining the problem domain according to the geometric dimensions of the plate and regularly arranging an appropriate number of nodes within the problem domain; The radial base point interpolation method is used to calculate the shape function of each node and approximate the displacement field u h (x) is expressed as: u h (x)={r p}G -1 U=ΦU Φ={rp}G -1 Among them, Φ is the shape function, r and p are radial basis functions and polynomial basis functions, and G is the basis function coefficient.

Images