A numerical solution for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations

Abstract

This invention discloses a numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations. The method involves constructing a minimum rectangular domain that can enclose the two-dimensional medium-thick plate, and transforming this domain into a standard square domain using a Jacobian matrix. Gaussian integration points are generated on the standard square domain. The geometric boundary of the two-dimensional medium-thick plate is defined. The elastic stiffness matrix of the two-dimensional medium-thick plate is calculated. A two-dimensional elastic stiffness matrix with the same dimension as the Gaussian integration points is formed. The stiffness at the Gaussian integration points within the open domain of the minimum rectangular domain is set to 0, forming a variable stiffness matrix. The energy functional variational problem of buckling of two-dimensional medium-thick plates is solved based on the energy principle and the Ritz method. The energy functional formula and solution procedure of the "global trial function + variable stiffness matrix" model proposed in this invention are completely standard, and the solution approach is applicable to other two-dimensional solid mechanics variational problems. The entire solution process is meshless and nodeless, freeing researchers from tedious finite element modeling and enabling more efficient solutions to engineering problems.

Description

Technical Field

[0001] This invention relates to the field of computational structural mechanics, and in particular to a numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations. Background Technology

[0002] Two-dimensional medium-thick plate structures are fundamental structural components widely used in engineering. Instability is one of the most typical failure modes leading to structural failure and damage. Ensuring the stability of two-dimensional medium-thick plates is crucial for reducing the risk of structural failure and ensuring structural safety, thus preventing overall structural instability. In engineering projects, the geometric configuration of two-dimensional medium-thick plates changes with design requirements under different working conditions. With these changes, stability analysis based on finite element method (FEM) technology requires remodeling according to the new geometric configuration. This often places a heavy workload on engineers. Furthermore, FEM numerical solution techniques based on mesh models are difficult to model and computationally intensive when dealing with intricate and complex geometric boundaries. Therefore, it is necessary to develop new meshless numerical solution methods for high-precision stability analysis of two-dimensional medium-thick plates with arbitrary geometric configurations. This can avoid the need for finite element mesh model reconstruction due to changes in the geometric model and improve the efficiency of numerical simulation. Summary of the Invention

[0003] To address the aforementioned technical problems, this invention proposes a numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations. This method can solve the buckling problem of two-dimensional medium-thick plates with arbitrary geometric configurations, and the solution process is meshless and nodeless, directly solving based on the geometric model, thus freeing engineers from tedious finite element modeling.

[0004] To achieve the above objectives, the technical solution of the present invention is as follows:

[0005] A numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations includes the following steps:

[0006] Obtain the material constants, thickness, geometric boundary equations, boundary condition parameters, and loads of a two-dimensional medium-thick plate with arbitrary geometric configuration; construct the minimum rectangular domain that can enclose the two-dimensional medium-thick plate, and transform the minimum rectangular domain into a standard square domain using the Jacobian matrix;

[0007] The two-dimensional medium-thick plate is discretized using Gaussian integral points in the interval [-1, 1] within a standard square domain;

[0008] Define the geometric boundary of a two-dimensional medium-thick plate according to the solid geometric boundary equation;

[0009] Based on the material constants of a two-dimensional medium-thick plate, the elastic stiffness matrix of the two-dimensional medium-thick plate is calculated. The elastic stiffness matrix includes an in-plane stiffness matrix A, a coupled stiffness matrix B, a bending stiffness matrix D, and a shear stiffness matrix A₀.s ;

[0010] Based on the Gaussian integration points and the elastic stiffness matrix, a two-dimensional elastic stiffness matrix with the same dimension as the Gaussian integration points is formed.

[0011] Based on the two-dimensional elastic stiffness matrix and the geometric boundary of the two-dimensional medium-thick plate, the stiffness at the Gaussian integral point in the opening domain of the minimum rectangular domain is set to 0, forming a variable stiffness matrix.

[0012] Based on the geometric boundary equations and boundary condition parameters, a global trial function considering the boundary conditions is constructed. Using the energy principle and Ritz method, pre-buckling analysis is performed on a two-dimensional medium-thick plate using the global trial function and the variable stiffness matrix within the minimum rectangular domain to obtain in-plane stress. Buckling analysis is then performed based on the in-plane stress, the global trial function, and the variable stiffness matrix within the minimum rectangular domain to obtain the buckling stiffness matrix eigenvalue equation. A Gaussian numerical integral is performed on each energy integral term in the buckling stiffness matrix eigenvalue equation, removing the energy in the opening region within the minimum rectangular domain. The eigenvalue problem is solved to obtain the buckling load and modes of the two-dimensional medium-thick plate.

[0013] Preferably, the in-plane stiffness matrix A, the coupled stiffness matrix B, the bending stiffness matrix D, and the shear stiffness matrix A... s As shown below:

[0014]

[0015] Each element is calculated according to the following formula:

[0016]

[0017] in denoted as the stiffness coefficient in the thickness direction, h is the thickness of the two-dimensional medium-thick plate, and z is the coordinate in the thickness direction.

[0018] Preferably, a two-dimensional elastic stiffness matrix with the same dimension as the Gaussian integration point is formed based on the Gaussian integration point and the elastic stiffness matrix, specifically including the following steps:

[0019] Based on the number of Gaussian integration points M and N in the directions ξ and η, the stiffness matrices A, B, D, and A are... s Each element A in ij B ij D ij Extended to a two-dimensional elastic stiffness matrix of dimension M×N:

[0020]

[0021] A total of 21 two-dimensional elastic stiffness matrices are generated, and all element values ​​in each of these 21 two-dimensional elastic stiffness matrices are the same.

[0022] Preferably, the step of setting the stiffness at the Gaussian integration point within the opening domain of the minimum rectangular domain to 0, based on the two-dimensional elastic stiffness matrix and the geometric boundary of the two-dimensional medium-thick plate, to form a variable stiffness matrix, specifically includes the following steps:

[0023] Based on the generated M×N Gaussian integration points, and combined with the defined geometric boundary of the two-dimensional medium-thick plate: determine whether each Gaussian integration point is located at the outer boundary Γ0 and the inner boundary Γ0. i Whether it is inside or outside, if the Gaussian integration point is located outside the outer boundary Γ0, or inside the inner boundary Γ0. i Within, that is, within the following set Inside:

[0024]

[0025] Then the point is determined to belong to the standard square domain Ω. s Within the opening, where e is the total number of inner boundaries of the two-dimensional medium-thick plate; according to equations (8) to (10), 21 two-dimensional elastic stiffness matrices are placed in the opening domain. The stiffness values ​​at the Gaussian integration points within the inner region are set to 0:

[0026]

[0027] The above equation is the variable stiffness matrix expression for simulating the geometric configuration and energy distribution of any two-dimensional medium-thick plate.

[0028] Preferably, based on the energy principle and the Ritz method, pre-buckling analysis is performed on a two-dimensional medium-thick plate using a global trial function and a variable stiffness matrix within a minimum rectangular domain to obtain in-plane stress. Then, buckling analysis is performed based on the in-plane stress, the global trial function, and the variable stiffness matrix within the minimum rectangular domain to obtain the eigenvalue equation of the buckling stiffness matrix. Specifically, the steps include:

[0029] The total potential energy during the pre-buckling stage is:

[0030] Π1=U m -U F (twenty two)

[0031] Among them U m The membrane strain energy is defined as:

[0032]

[0033] Where, ε 0 Let |J| be the linear strain in the mid-surface, |J| be the Jacobian determinant, and |J| = ab / 4.

[0034] U F The potential energy of the external load is defined as:

[0035]

[0036] Where a and b are the length and width of the minimum rectangular domain, respectively, and N(ξ,η) is the in-plane distributed load applied to the boundary Γ0, which needs to be decomposed into the directions ξ and η, i.e., N ξ =N(ξ,η)cosα(ξ,η) and N η =N(ξ,η)cosα(ξ,η),

[0037] Substituting the first two terms u0 and v0 of the global trial function into formulas (24) and (23), and then substituting the result into formula (22), we can find the functional Π1 with respect to the unknown coefficients. Resident value:

[0038]

[0039] This will give us a system of linear equations with unknown coefficients. By writing them in matrix form and inverting them, we can obtain the unknown coefficients:

[0040]

[0041] The vector W is derived from the load potential energy U. F Decision, K m The in-plane stiffness matrix is ​​shown below:

[0042]

[0043] unknown coefficients Substituting the global trial function into the equation, we obtain the in-plane displacements u0 and v0. Then, we obtain the mid-plane linear strain ε using equation (14). 0 , and then ε 0 Substituting the constitutive equation considering only in-plane stress, we get the following:

[0044]

[0045] Buckling analysis: As the external load increases, the out-of-plane deformation state caused by the aforementioned in-plane initial stress is called buckling. The total potential energy of buckling is as follows:

[0046] Π2=U b -U N (29)

[0047] Among them U N The strain energy caused by the in-plane initial stress is shown below:

[0048]

[0049] in For von Kármán nonlinear strain,

[0050] Out-of-plane bending strain energy U bfor:

[0051]

[0052] Where κ and γ are the curvature and transverse shear strain of the plate, respectively.

[0053] By finding the unknown coefficients of the functional Π2 The stationary value is:

[0054]

[0055] The eigenvalue equations for the buckling stiffness matrix are obtained as follows:

[0056] {K b -λK g}{C}=0 (33)

[0057] in The vector represents unknown coefficients, and each component contains J×K elements; the bending stiffness matrix K b and geometric matrix K g They are respectively:

[0058]

[0059]

[0060] They are all (3×J×K)×(3×J×K) matrices.

[0061] In formula (27), the in-plane stiffness matrix K m The bending stiffness matrix K of formula (34) b In this context, the matrix elements are defined as follows:

[0062]

[0063]

[0064] Where the matrix and Each element in and ( d, e, f, g = 0, 1) are defined as follows:

[0065]

[0066] and

[0067]

[0068] in χ j(ξ) and χ k (η) represents the j-th and k-th terms of the Legendre polynomial. and For Legendre polynomials, the first Item and the Item, element and Located in the matrix respectively and The r-th row and s-th column,

[0069] The geometric matrix K of formula (35) g The matrix elements are defined as follows:

[0070]

[0071] Where the matrix (i,j = 1,2,6; p,q = 0,1); Elements in d,e,f,g=0,1) Defined as:

[0072]

[0073] Where w r (ξ,η) and w s (ξ,η) is set in formula (39) It can be obtained, and Defined as:

[0074]

[0075] From the above formulas (36), (37), and (40), it can be seen that the variable stiffness matrix A of a two-dimensional medium-thick plate with arbitrary geometric configuration is... ij B ij D ij To pass through the standard square domain Ω s The inner integral characterizes the energy in the open region of the minimum rectangular domain, because the variable stiffness matrix A in equation (12) ij B ij D ij The stiffness at the Gaussian integration point within the opening region of the smallest rectangular domain has been set to 0.

[0076] Based on the above technical solution, the beneficial effects of this invention are as follows: This invention proposes a numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations. It simulates two-dimensional medium-thick plates with arbitrary geometric configurations by drilling holes in a standard geometric domain. Then, it constructs an elastic variable stiffness matrix to characterize the material distribution in space and discretizes it using Gaussian integral points. The energy functional is solved using global trial functions, thereby simulating the deformation of two-dimensional medium-thick plates with arbitrary geometric configurations. This solves the problem of difficulty in solving variational problems on complex geometric domains using global trial functions. The energy functional formula and solution procedure of the "global trial function + variable stiffness matrix" model proposed in this invention are completely standard, and the solution approach is applicable to other two-dimensional solid mechanics variational problems. The entire solution process is meshless and nodeless, freeing researchers from tedious finite element modeling and enabling more efficient solutions to engineering problems. Attached Figure Description

[0077] Figure 1 This is a definition of a two-dimensional medium-thick plate geometric model, load, discrete model, and ply definition, where: (a) a two-dimensional medium-thick plate geometric model and discrete Ritz method modeling; (b) two-dimensional medium-thick plate load decomposition; (c) two-dimensional medium-thick plate discrete model, where gray Gaussian integral points are located within the real plate domain Ω, and black Gaussian integral points are located within the opening domain with stiffness set to 0; (d) two-dimensional medium-thick plate ply and ply orientation angle definition.

[0078] Figure 2 The variable stiffness matrix representation and discretization of a two-dimensional medium-thick plate is based on Gaussian integral points, where: (a) a perforated medium-thick plate; (b) a variable stiffness matrix representation and discretization of a perforated medium-thick plate.

[0079] Figure 3 This is a flowchart of a numerical solution method for buckling of a two-dimensional medium-thick plate with arbitrary geometric configuration in one embodiment;

[0080] Figure 4 The buckling of a circular plate under distributed load, wherein: (a) the geometric model of the circular plate and the load; (b) the discrete model of the circular plate;

[0081] Figure 5 These are the first five buckling load coefficients and modal diagrams of a circular plate under simply supported and fixed boundary conditions, subjected to different distributed loads. Detailed Implementation

[0082] Based on the Ritz method, a numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations is proposed. This method can solve the stability problem of two-dimensional medium-thick plates with arbitrary geometric configurations, and the solution process is meshless and node-free, freeing engineers from tedious finite element modeling. The specific steps of this method are as follows:

[0083] Step 1: Construct a standard two-dimensional geometric solution domain. Input material constants, the thickness of the two-dimensional medium-thick plate, the structural geometric boundary equations, boundary conditions, and loads. Material constants include E1, E2, and G. 12 ν 12 The thickness of a two-dimensional medium-thick plate is h, such as... Figure 1 As shown in (a), the light gray region Ω represents the true geometric domain of a two-dimensional medium-thick plate. First, a minimum rectangular domain Ω that can cover the true geometric domain Ω of the two-dimensional medium-thick plate is constructed. r The length and width of the minimum rectangular area are a and b, respectively, and are determined by the following formula:

[0084] ξ=2x / a, η=2y / bξ, η∈[-1,1](1)

[0085] The smallest rectangular area Ω r Transform to standard square domain Ω s In the middle, Ω s ={(ξ,η)|-1≤ξ≤1,-1≤η≤1}, where ξ and η are the normalized coordinates in the x and y directions, respectively. This step transforms the minimum rectangular solution domain covering a two-dimensional medium-thick plate with arbitrary geometric configuration into a standard square solution domain.

[0086] Step 2: Generate Gaussian integration points on the standard square domain. The standard square domain Ω is generated using the Gauss-Legendary quadrature formula. s The Gaussian integration points and their weight coefficients within the standard square domain are determined. M and N Gaussian integration points are generated in the ξ and η directions, respectively, and then discretized.

[0087] Step 3: Define the inner and outer boundaries of the true geometric domain Ω of the two-dimensional medium-thick plate. Based on the input geometric boundary equations of the two-dimensional medium-thick plate, denote the outer boundary of the design domain Ω as Γ0 and the inner boundary as Γ0. i Where i represents the i-th inner boundary, such as Figure 1 As shown in (a), the inner and outer boundaries are functions of the coordinates (ξ, η), and their boundary function equations can be expressed using level set functions as follows:

[0088] Γ0(ξ,η)=0,Γ i (ξ, η)=0, i=1,2,3,... (2)

[0089] Therefore, whether a point (ξ,η) is inside or outside the boundary can be determined by the following formula:

[0090]

[0091] Step 4: Calculate the elastic stiffness matrix of the two-dimensional medium-thick plate. Here, the elastic stiffness of anisotropic materials is uniformly used for calculation. Based on the input material constants, the elastic stiffness matrix of the k-th layer in the thickness direction of the two-dimensional medium-thick plate is calculated. k Elastic stiffness:

[0092]

[0093]

[0094] in Let be the transfer reduced stiffness coefficient of the k-th layer, when the material is isotropic. The same in each layer is θ k For the pavement orientation of the k-th layer, see... Figure 1 (d); U i (i = 1, 2, ..., 5) represents the material invariants:

[0095]

[0096] Q ij (i,j=1,2,6;4,5) are the elastic constants of the composite material:

[0097]

[0098] Among them, E1, E2, G 12 ν 12 Let ν be a material constant, and κ be a shear correction factor κ = 5 / 6. When the material is isotropic, the elastic modulus is E = E1 = E2, and the Poisson's ratio is ν = ν 12 =ν 21 And the shear modulus G = G 12 =G 13 =G 23 =E / [2(1+ν)].

[0099] Calculate the in-plane stiffness matrix A, coupled stiffness matrix B, bending stiffness matrix D, and shear stiffness matrix A according to formulas (4)–(7). s :

[0100]

[0101] Each element is calculated according to the following formula:

[0102]

[0103] in denoted as the stiffness coefficient in the thickness direction, h is the thickness of the two-dimensional medium-thick plate, and z is the coordinate in the thickness direction.

[0104] Step 5: Generate a two-dimensional elastic stiffness matrix based on Gaussian integration points and elements of the stiffness matrix. According to the number of Gaussian integration points M and N in the ξ and η directions, respectively, the stiffness matrices A, B, D, and A... s Each element A in ij Bij D ij Extended to a two-dimensional elastic stiffness matrix of dimension M×N:

[0105]

[0106] The above operations generate a total of 21 two-dimensional elastic stiffness matrices, and all element values ​​in each of these 21 two-dimensional elastic stiffness matrices are identical.

[0107] Step 6: Generate a variable stiffness matrix based on the two-dimensional elastic stiffness matrix and the two-dimensional medium-thick plate geometric boundary. Based on the M×N Gaussian integration points generated in Step 2, and combining equation (3): For each Gaussian integration point, first determine whether the Gaussian integration point is located at the outer boundary Γ0 and the inner boundary Γ0. i Is it inside or outside? If the Gaussian integration point is located outside the outer boundary Γ0, or inside the inner boundary Γ0... i Within, that is, within the following set Inside:

[0108]

[0109] Then it is determined that the point belongs to the opening of the standard square domain, where e is the total number of the inner boundaries of the two-dimensional medium-thick plate. According to equations (8) to (10), the 21 elastic stiffness matrices are placed in the opening domain. The stiffness values ​​at the Gaussian integration points within the inner region are set to 0:

[0110]

[0111] The above equation is the variable stiffness matrix expression for simulating the geometric configuration and energy distribution of an arbitrary two-dimensional medium-thick plate. For example... Figure 1 As shown in (c), the black dots represent Gaussian integral points with stiffness of 0, located in the open-hole domain. Inside; gray dots represent Gaussian integral points with stiffness, located within the geometric domain Ω of a real two-dimensional medium-thick plate.

[0112] To illustrate the variable stiffness matrix of a two-dimensional medium-thick plate, for Figure 2 (a) shows a medium-thick plate with an opening. The number of Gaussian integration points in the ξ and η directions are M = N = 10, respectively. First, the stiffness matrix is ​​expanded to a 10×10 two-dimensional elastic stiffness matrix. Then, based on the outer boundary Γ0 and the inner boundary Γ... i The horizontal set function sets the stiffness value within the opening to 0, such as... Figure 2 As shown in (b), the variable stiffness matrix representation and discretization of a two-dimensional medium-thick plate with arbitrary geometric configuration are completed. Therefore, for the 21 stiffness values ​​in equation (10), a total of 21 variable stiffness matrices A are generated. ij B ij D ij(i,j=1,2,6;4,5), as shown in Equation (12). These variable stiffness matrices can be used to characterize the strain energy of a two-dimensional medium-thick plate with arbitrary geometric configuration by multiplying the Hadamard product with the strain matrix during subsequent numerical integration.

[0113] The above six steps can be used to characterize and discretize the variable stiffness of a two-dimensional medium-thick plate with arbitrary geometric configuration.

[0114] Step 7: Solve the functional variational problem of buckling energy in two-dimensional medium-thick plates based on the Ritz method and first-order shear deformation theory. This includes the following steps:

[0115] Step 701: Based on the first-order shear deformation theory and Legendre polynomials, construct a global trial function considering boundary conditions. According to the first-order shear deformation theory, in the Cartesian coordinate system, assume that the deformations u, v, w at any point in the medium-thick plate along the x, y, and z directions are the mid-surface displacements u0, v0, w0 and the rotation angle ψ. u and ψ v Functions:

[0116]

[0117] Where ψ u and ψ v These are the rotation angles of the sections about the y-axis and x-axis, respectively.

[0118] The strain of a medium-thick plate is defined as:

[0119]

[0120]

[0121] Where, ε 0 ,κ, γ represents the mid-surface linear strain, curvature, von Kármán nonlinear strain, and transverse shear strain, respectively. The unknowns u0, v0, w0, and rotation angle ψ are solved. u and ψ v It can draw the buckling deformation modes of medium-thick plates.

[0122] Construct global trial functions u0, v0, w0, ψ considering boundary conditions. x ψ y :

[0123]

[0124] in The coefficients are unknown; J and K represent the number of terms in the Legendre polynomials in the directions ξ and η, respectively; χ j (ξ), χ k (η) are the Legendre polynomials in the directions ξ and η, respectively:

[0125]

[0126] Where (j-1) is the exponent of the Legendre polynomial. They are the boundary Γ i Displacement functions u, v, w, ψ on (i = 0, 1, 2, ..., e) u ,ψ v The boundary condition function. If Γ i The p-th edge of (i = 0, 1, 2, ..., e) is Z ip (ξ,η), then it can be obtained through the boundary function The edge Z ip Apply an exponent to (ξ,η) Characterizing boundary condition functions (θ=u,v,w,ψ u ,ψ v ):

[0127]

[0128] Where V i For the boundary Γ i The total number of edges.

[0129] For example, the boundary conditions for a rectangular plate can be set at the outer boundary. Boundary indices are applied to the four edges in a counter-clockwise direction:

[0130]

[0131] For Γ i A displacement component of the p-th edge of (i = 0, 1, 2, ..., e) Loading constraints, boundary indices It can be defined by the following formula:

[0132]

[0133] For free, simply supported, and fixed boundary conditions, the exponents of the five displacement components on the p-th edge are... The definition is as follows:

[0134]

[0135] Therefore, by multiplying the above four formulas and multiple boundary condition equations, we obtain the boundary condition equations for a two-dimensional medium-thick plate.

[0136] Step 702: Pre-buckling and buckling analysis to obtain in-plane stress distribution and buckling eigenvalue equations. Based on the energy principle and Ritz method, pre-buckling analysis is performed on a two-dimensional medium-thick plate using a global trial function and a variable stiffness matrix within a minimum rectangular domain to obtain in-plane stress. Buckling analysis is then performed based on the in-plane stress, the global trial function, and the variable stiffness matrix within the minimum rectangular domain to obtain the buckling stiffness matrix eigenvalue equations. Taking the buckling of a two-dimensional medium-thick plate as an example, the total potential energy in the pre-buckling stage is:

[0137] Π1=U m -U F (twenty two)

[0138] Among them U m The membrane strain energy is defined as:

[0139]

[0140] Where, ε 0 Let ε be the linear strain in the mid-surface, |J| be the Jacobian determinant and |J| = ab / 4, where ε 0 For the mid-surface linear strain.

[0141] U F The potential energy of the external load is defined as:

[0142]

[0143] Where N(ξ,η) is the in-plane distributed load applied to the boundary Γ0, and it needs to be decomposed into the directions ξ and η, i.e., N ξ =N(ξ,η)cosα(ξ,η) and N η =N(ξ,η)cosα(ξ,η), such as Figure 1 As shown in (b).

[0144] Substitute the first two terms of formula (16) into formulas (24) and (23), and then substitute the result into formula (22). Next, calculate the functional Π1 with respect to the unknown coefficients. Resident value:

[0145]

[0146] This will give us a system of linear equations with unknown coefficients. By writing them in matrix form and inverting them, we can obtain the unknown coefficients:

[0147]

[0148] The vector W is derived from the load potential energy U. F Decision, K m Here is the in-plane stiffness matrix:

[0149]

[0150] unknown coefficients Substituting into formula (16) yields the in-plane displacements u0 and v0, and then the mid-plane linear strain ε is obtained through formula (14). 0 , and then ε 0 Substituting the constitutive equation considering only in-plane forces:

[0151]

[0152] Buckling analysis is performed below. As the external load increases, the out-of-plane deformation state caused by the aforementioned in-plane initial stress is called buckling. The total potential energy of buckling is:

[0153] Π2=U b -U N (29)

[0154] Among them U N The strain energy caused by the in-plane initial stress:

[0155]

[0156] in For von Kármán nonlinear strain.

[0157] The out-of-plane bending strain energy U b for:

[0158]

[0159] Where κ and γ are the curvature and transverse shear strain of the plate, respectively.

[0160] By finding the unknown coefficients of the functional Π2 Resident value:

[0161]

[0162] The eigenvalue equations for the buckling rigid matrix are obtained as follows:

[0163] {K b -λK g}{C}=0 (33)

[0164] in The vector represents unknown coefficients, and each component contains J×K elements; the bending stiffness matrix K b and geometric matrix K g They are respectively:

[0165]

[0166]

[0167] They are all (3×J×K)×(3×J×K) matrices.

[0168] In formula (27), the in-plane stiffness matrix K m The bending stiffness matrix K of formula (34) b In this context, the matrix elements are defined as follows:

[0169]

[0170]

[0171] Where the matrix and Each element in and ( d, e, f, g = 0, 1) are defined as follows:

[0172]

[0173] and

[0174]

[0175] in χ j (ξ) and χ k (η) represents the j-th and k-th terms of the Legendre polynomial. and For Legendre polynomials, the first Item and the Item. Element and Located in the matrix respectively and The r-th row and s-th column.

[0176] The geometric matrix K of formula (35) g The matrix elements are defined as follows:

[0177]

[0178] Where the matrix (i,j = 1,2,6; p,q = 0,1); Elements in d,e,f,g=0,1) Defined as:

[0179]

[0180] Where w r (ξ,η) and w s (ξ,η) is set in formula (39) It can be obtained, and Defined as:

[0181]

[0182] From the above formulas (36), (37), and (40), it can be seen that the variable stiffness matrix A of a two-dimensional medium-thick plate with arbitrary geometric configuration is... ij B ij D ij To pass through the standard square domain Ω s The inner integral characterizes the energy in the open region of the minimum rectangular domain, because the variable stiffness matrix A in equation (12) ij B ij D ij The stiffness at the Gaussian integration point within the opening region of the smallest rectangular domain has been set to 0.

[0183] Step 703: Perform a Gaussian numerical integration based on the variable stiffness matrix for each energy integral term in the buckling stiffness matrix eigenvalue equation to remove the energy in the opening region. When equations (38) and (41) are incorporated into the double integral and Gaussian numerical integration is performed: Perform a Gaussian numerical integration for each term in the energy integral to remove the energy in the opening region:

[0184]

[0185] Energy has been removed from the open region, where g(ξ,η) represents the integrand in equations (38) and (41), and ξ m η n These represent the Gaussian integration points in the interval [-1, 1] along the directions ξ and η, respectively, and their weight coefficients are respectively and In the above formula, the number of Gaussian integration points is usually large because a large number of Gaussian integration points are needed to capture the geometric boundaries of complex opening regions, so as to accurately calculate the energy of the opening regions.

[0186] The above three steps complete the numerical integration energy removal of the complex geometric domain. Finally, by solving the eigenvalue problem in formula (33), the buckling load coefficient and modes of a two-dimensional medium-thick plate with arbitrary geometric configuration can be obtained. The flowchart of the discrete energy method for two-dimensional medium-thick plates is summarized in [link to flowchart]. Figure 3 middle.

[0187] Example

[0188] To overcome the geometric contradictions of existing numerical methods, a numerical solution method for buckling two-dimensional medium-thick plates with arbitrary geometric configurations—the discrete energy method—is proposed. This embodiment analyzes the buckling of a circular plate under different boundary conditions and distributed loads. The distributed load angle is defined as β, as shown below. Figure 4 As shown in (a), the radius of the circular plate is R = 5m and the thickness is h = 0.1m.

[0189] The discrete energy method is applied to the buckling problem of a circular plate under distributed loads, considering both simply supported and fixed boundary conditions, to solve for the buckling load coefficient and modes of the circular plate. The material constant of the circular plate is: E = 2.0 × 10⁻⁶. 11 N / m 2 Poisson's ratio υ = 0.3, density ρ = 8000 kg / m³ 3 Bending stiffness is defined as D = Eh 3 / [12(1-ν 2 The outer boundary Γ0 of the circular plate is defined as:

[0190]

[0191] Where a and b are the length and width of the smallest rectangular region covering the circular plate, respectively. The set of opening regions within the smallest rectangular region. Defined as:

[0192]

[0193] For both simply supported and fixed-supported boundary conditions, the boundary conditions of the outer boundary... Defined as:

[0194]

[0195] This embodiment only considers two boundary conditions: simply supported and fixed.

[0196] The geometric model of the circular plate is as follows Figure 4 As shown: Figure 4 In (a) and (b), the minimum rectangular domain covering the circular plate has dimensions a = b = 10 m, an inner hole radius R = 5 m, and a plate thickness h = 0.1 m. The buckling load coefficient and modes of the circular plate under simply supported and fixed boundary conditions were analyzed using the discrete energy method. The structure was discretized in the ξ and η directions using M = N = 500 Gaussian integration points, as shown below. Figure 4 As shown in (b), the opening area The Gaussian integration points within the circular plate domain are black, and those within the circular plate domain are gray. The number of terms in the ξ and η directional shape functions are J = K = 10 terms each. The specific implementation steps are as follows:

[0197] Step 1: Construct a standard two-dimensional geometric solution domain. Input material constants, the thickness of the two-dimensional medium-thick plate, the structural geometric boundary equations, boundary conditions, and loads. Input material constants E1, E2, and G. 12 ν 12 The thickness h of a two-dimensional medium-thick plate. For example... Figure 4 As shown in (a), the light gray region Ω represents the true geometric domain of a two-dimensional medium-thick plate. First, a minimum rectangular domain Ω that can cover the true geometric domain Ω of the two-dimensional medium-thick plate is constructed.r The length and width of the minimum rectangular area are a and b, respectively, and are determined by the following formula:

[0198] ξ=2x / a, η=2y / bξ, η∈[-1,1](48)

[0199] The smallest rectangular area Ω r Transform to standard square domain Ω s In the middle, Ω s ={(ξ,η)|-1≤ξ≤1,-1≤η≤1}, where ξ and η are the normalized coordinates in the x and y directions, respectively. This step transforms the minimum rectangular solution domain covering a two-dimensional medium-thick plate with arbitrary geometric configuration into a standard square solution domain Ω. s middle.

[0200] Step 2: Generate Gaussian integration points on the standard square domain. The standard square domain Ω is generated using the Gauss-Legendary quadrature formula. s The Gaussian integration points and their weight coefficients within the standard square domain are generated. M×N = 500×500 Gaussian integration points are generated in the ξ and η directions within the domain and discretized, as shown below. Figure 4 As shown in (b).

[0201] Step 3: Define the outer boundary of the true geometric domain Ω of the two-dimensional medium-thick circular plate. Based on the input geometric boundary equation of the two-dimensional medium-thick plate, denote the outer boundary of the design domain Ω as Γ0, as follows: Figure 4 As shown in (a), the outer boundary is a function of the coordinates (ξ, η), and its boundary function equation can be expressed using the level set function as follows:

[0202]

[0203] Therefore, whether a point (ξ,η) is inside or outside the boundary can be determined by the following formula:

[0204]

[0205] According to the above equation, in Figure 4 In (b), the black dots are Gaussian integration points with stiffness set to 0, and the gray dots are Gaussian integration points with non-zero stiffness.

[0206] Step 4: Calculate the elastic stiffness matrix of the two-dimensional medium-thick plate. For the isotropic material used in this embodiment, the elastic stiffness of anisotropic materials is uniformly used for calculation. Therefore, the material constants of isotropic materials are defined as: elastic modulus E1 = E2 = E, Poisson's ratio ν. 12 =ν 21 =ν, shear modulus G 12 =G=E / [2(1+ν)]. According to formulas (4) to (8) in the steps, we can calculate:

[0207]

[0208] Step 5: Generate a two-dimensional elastic stiffness matrix based on Gaussian integration points and the elastic stiffness matrix. According to the number of Gaussian integration points in the ξ and η directions, the above stiffness matrices A, B, D, and A... s Each element A in ij B ij D ij Expanded to a two-dimensional elastic stiffness matrix with dimensions 500×500:

[0209]

[0210] The above operations generate a total of 21 two-dimensional elastic stiffness matrices, and all element values ​​in each of these 21 two-dimensional elastic stiffness matrices are identical.

[0211] Step 6: Generate a variable stiffness matrix based on the two-dimensional elastic stiffness matrix and the two-dimensional medium-thick plate geometric boundary. Based on the 500×500 Gaussian integration points generated in Step 2, combined with Equation (49): For each Gaussian integration point, first determine whether the Gaussian integration point is located inside or outside the outer boundary Γ0. If the Gaussian integration point is located outside the outer boundary Γ0, that is, within the following set... Inside:

[0212]

[0213] Then it is determined that the point belongs to the opening of the standard square domain. According to equations (49) to (53), the 21 two-dimensional elastic stiffness matrices are placed in the opening domain. The stiffness values ​​at the Gaussian integration points within the inner region are set to 0:

[0214]

[0215] The above equation is the variable stiffness matrix expression for simulating the geometric configuration and energy distribution of an arbitrary two-dimensional medium-thick plate. For example... Figure 4 As shown in (b), the black Gaussian integration points represent Gaussian integration points with stiffness of 0, located in the open-hole domain. Inside; while the gray Gaussian integral points represent Gaussian integral points with stiffness, located within the geometric domain Ω of a real two-dimensional medium-thick circular plate.

[0216] Step 7: Solve the functional variational problem of buckling energy in two-dimensional medium-thick plates based on the Ritz method and first-order shear deformation theory. Based on the variable stiffness matrix A constructed above... ij B ij D ijTo solve the functional stationary value problem, see formulas (13) to (43). The variable stiffness matrix (54) needs to be substituted into formulas (36) to (43) to remove the energy of the opening region in the minimum rectangular domain. Then the buckling load coefficient and mode of the circular plate under distributed load can be obtained.

[0217] The buckling load coefficients and modes of a circular plate under simply supported and fixed boundary conditions were solved using the proposed discrete energy method. The results are shown in Table 1 and compared with existing results. It can be seen that the results of the discrete energy method analysis are in good agreement with existing results. Figure 5 The first five buckling modes of a circular plate under two different boundary conditions and different distributed loads are presented.

[0218] Table 1 Dimensionless buckling load coefficient λ of circular plate under loads distributed at different angles β cr =N cr R 2 / D

[0219]

[0220]

[0221] It should be understood that although the steps in the flowchart above are shown sequentially as indicated by the arrows, these steps are not necessarily executed in the order indicated by the arrows. Unless explicitly stated herein, there is no strict order restriction on the execution of these steps, and they can be executed in other orders. Moreover, at least some steps in the flowchart above may include multiple sub-steps or multiple stages. These sub-steps or stages are not necessarily completed at the same time, but can be executed at different times. The execution order of these sub-steps or stages is not necessarily sequential, but can be performed alternately or in turn with other steps or at least some of the sub-steps or stages of other steps.

[0222] The above description is merely a preferred embodiment of the numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations disclosed in this invention, and is not intended to limit the scope of protection of the embodiments in this specification. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the embodiments in this specification should be included within the scope of protection of the embodiments in this specification.

Claims

1. A numerical solution method for buckling of two-dimensional medium-thick plates with arbitrary geometric configurations, characterized in that, Includes the following steps: Obtain the material constants, thickness, geometric boundary equations, boundary condition parameters, and loads of a two-dimensional medium-thick plate with arbitrary geometric configuration; construct the minimum rectangular domain that can enclose the two-dimensional medium-thick plate, and transform the minimum rectangular domain into a standard square domain using the Jacobian matrix; The two-dimensional medium-thick plate is discretized using Gaussian integral points in the interval [-1, 1] within a standard square domain; The geometric boundary of a two-dimensional medium-thick plate is defined based on the geometric boundary equation of the two-dimensional medium-thick plate. Based on the material constants of a two-dimensional medium-thick plate, the elastic stiffness matrix of the two-dimensional medium-thick plate is calculated. The elastic stiffness matrix includes an in-plane stiffness matrix A, a coupled stiffness matrix B, a bending stiffness matrix D, and a shear stiffness matrix A₀. s ; Based on the Gaussian integration points and the elastic stiffness matrix, a two-dimensional elastic stiffness matrix with the same dimension as the Gaussian integration points is formed. Based on the two-dimensional elastic stiffness matrix and the geometric boundary of the two-dimensional medium-thick plate, the stiffness at the Gaussian integral point in the opening domain of the minimum rectangular domain is set to 0, forming a variable stiffness matrix. Construct a global trial function considering the boundary conditions based on the parameters of the geometric boundary equations and boundary conditions; Based on the energy principle and the Ritz method, pre-buckling analysis is performed on a two-dimensional medium-thick plate using a global trial function and a variable stiffness matrix within a minimum rectangular domain to obtain in-plane stress. Buckling analysis is then performed based on the in-plane stress, the global trial function, and the variable stiffness matrix within the minimum rectangular domain to obtain the eigenvalue equation of the buckling stiffness matrix. For each energy integral term in the eigenvalue equation of the buckling stiffness matrix, a Gaussian numerical integral is performed based on the variable stiffness matrix, removing the energy in the opening region within the minimum rectangular domain. Solving the eigenvalue problem yields the buckling load and modes of the two-dimensional medium-thick plate.

2. The numerical solution method for buckling of a two-dimensional medium-thick plate with arbitrary geometric configuration according to claim 1, characterized in that, The in-plane stiffness matrix A, coupled stiffness matrix B, bending stiffness matrix D, and shear stiffness matrix A s As shown below: Each element is calculated according to the following formula: in denoted as the stiffness coefficient in the thickness direction, h is the thickness of the two-dimensional medium-thick plate, and z is the coordinate in the thickness direction.

3. The numerical solution method for buckling of a two-dimensional medium-thick plate with arbitrary geometric configuration according to claim 1, characterized in that, The step of forming a two-dimensional elastic stiffness matrix with the same dimension as the Gaussian integral points based on the Gaussian integral points specifically includes the following steps: Based on the number of Gaussian integration points M and N in the directions ξ and η, the stiffness matrices A, B, D, and A are... s Each element A in ij B ij D ij Extended to a two-dimensional elastic stiffness matrix of dimension M×N: A total of 21 two-dimensional elastic stiffness matrices are generated, and all element values ​​in each of these 21 two-dimensional elastic stiffness matrices are the same.

4. The numerical solution method for buckling of a two-dimensional medium-thick plate with arbitrary geometric configuration according to claim 1, characterized in that, The method based on the two-dimensional elastic stiffness matrix and the geometric boundary of the two-dimensional medium-thick plate, sets the stiffness at the Gaussian integration point within the opening domain of the minimum rectangular domain to 0, forming a variable stiffness matrix, specifically includes the following steps: Based on the generated M×N Gaussian integration points, and combined with the defined geometric boundary of the two-dimensional medium-thick plate: determine whether each Gaussian integration point is located at the outer boundary Γ0 and the inner boundary Γ0. i Whether it is inside or outside, if the Gaussian integration point is located outside the outer boundary Γ0, or inside the inner boundary Γ0. i Within, that is, within the following set Inside: Then the point is determined to belong to the standard square domain Ω. s Within the opening, where e is the total number of inner boundaries of the two-dimensional medium-thick plate; according to equations (8) to (10), 21 two-dimensional elastic stiffness matrices are placed in the opening domain. The stiffness values ​​at the Gaussian integration points within the inner region are set to 0: The above equation is the variable stiffness matrix expression for simulating the geometric configuration and energy distribution of any two-dimensional medium-thick plate.

5. The numerical solution method for buckling of a two-dimensional medium-thick plate with arbitrary geometric configuration according to claim 1, characterized in that, Based on the energy principle and Ritz method, pre-buckling analysis is performed on a two-dimensional medium-thick plate using a global trial function and a variable stiffness matrix within a minimum rectangular domain to obtain in-plane stress. Then, buckling analysis is performed based on the in-plane stress, the global trial function, and the variable stiffness matrix within the minimum rectangular domain to obtain the eigenvalue equation of the buckling stiffness matrix. The specific steps include: The total potential energy during the pre-buckling stage is: Π1=U m -U F (22) Among them U m The membrane strain energy is defined as: Where, ε 0 Let |J| be the linear strain in the mid-surface, |J| be the Jacobian determinant, and |J| = ab / 4. U F The potential energy of the external load is defined as: Where a and b are the length and width of the minimum rectangular domain, respectively, and N(ξ,η) is the in-plane distributed load applied to the boundary Γ0, which needs to be decomposed into the directions ξ and η, i.e., N ξ =N(ξ,η)cosα(ξ,η) and N η =N(ξ,η)cosα(ξ,η), Substituting the first two terms u0 and v0 of the global trial function into formulas (24) and (23), and then substituting the result into formula (22), we can find the functional Π1 with respect to the unknown coefficients. Resident value: This will give us a system of linear equations with unknown coefficients. By writing them in matrix form and inverting them, we can obtain the unknown coefficients: The vector W is derived from the load potential energy U. F Decision, K m The in-plane stiffness matrix is ​​shown below: unknown coefficients Substituting the global trial function into the equation, we obtain the in-plane displacements u0 and v0. Then, we obtain the mid-plane linear strain ε using equation (14). 0 , and then ε 0 Substituting the constitutive equation considering only in-plane stress, we get the following: Buckling analysis: As the external load increases, the out-of-plane deformation state caused by the aforementioned in-plane initial stress is called buckling. The total potential energy of buckling is as follows: Π2=U b -U N (29) Among them U N The strain energy caused by the in-plane initial stress is shown below: in For von Kármán nonlinear strain, Out-of-plane bending strain energy U b for: Where κ and γ are the curvature and transverse shear strain of the plate, respectively. By finding the unknown coefficients of the functional Π2 The stationary value is: The eigenvalue equations for the buckling stiffness matrix are obtained as follows: {K b -λK g }{C}=0 (33) in The vector represents unknown coefficients, and each component contains J×K elements; the bending stiffness matrix K b and geometric matrix K g They are respectively: They are all (3×J×K)×(3×J×K) matrices. In formula (27), the in-plane stiffness matrix K m The bending stiffness matrix K of formula (34) b In this context, the matrix elements are defined as follows: Where the matrix and Each element in and They are defined as follows: and in χ j (ξ) and χ k (η) represents the j-th and k-th terms of the Legendre polynomial. and For Legendre polynomials, the first Item and the Item, element and Located in the matrix respectively and The r-th row and s-th column, The geometric matrix K of formula (35) g The matrix elements are defined as follows: Where the matrix elements in Defined as: Where w r (ξ,η) and w s (ξ,η) is set in formula (39) It can be obtained, and Defined as: From the above formulas (36), (37), and (40), it can be seen that the variable stiffness matrix A of a two-dimensional medium-thick plate with arbitrary geometric configuration is... ij B ij D ij To pass through the standard square domain Ω s The inner integral characterizes the energy in the open region of the minimum rectangular domain, because the variable stiffness matrix A in equation (12) ij B ij D ij The stiffness at the Gaussian integration point within the opening region of the smallest rectangular domain has been set to 0.

Images