A method for quickly and stably solving a post-buckling path with multiple critical points based on a generalized displacement control method

Abstract

This invention discloses a method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method, comprising: S1, initializing state variables and anchored secant predictors; S2, predicting reference displacement, calibration displacement, and load increment factor based on state variables and anchored secant predictors; S3, calculating the displacement increment of the current iteration step based on the reference displacement, calibration displacement, and load increment factor, and updating the total structural displacement and total load; S4, updating the structural geometry and element internal forces based on the displacement increment of the current iteration step; S5, calculating the difference between the total structural internal forces and the external load based on the updated structural geometry and element internal forces to obtain the structural unbalanced force, and performing error verification based on the structural unbalanced force; S6, cyclically executing S1-S5 until the structural unbalanced force meets the equilibrium condition, and starting the next incremental step; the incremental analysis is executed until the user requirements are met, and the buckling path after multiple critical points is output.

Description

Technical Field

[0001] This invention belongs to the field of nonlinear finite element and structural nonlinear analysis technology, and particularly relates to a method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method. Background Technology

[0002] Over the past few decades, with advancements in materials preparation technology, engineers have focused on integrating high-performance materials with structures to achieve both lightweighting and improved structural strength and durability. However, when such thin-walled components are subjected to abnormal or extreme loads, they often experience significant nonlinear deformation due to structural instability. To meet the design needs of thin-walled structures in various engineering projects, including high-rise buildings, aerospace, machinery manufacturing, and bridge structures, many design codes require stability analysis as an integral part of structural design practices.

[0003] Traditionally, stability analysis has only been applied to simple structures; for complex structures, approximation techniques are mostly used. To ensure structural reliability, design codes typically require high safety factors, which significantly increases material consumption. Today, with the development of finite element method (FEM) technology, engineers can gain a more realistic understanding of the mechanical response of structures. As my country's infrastructure construction shifts from a phase of incremental development to high-quality development, the requirements for structural safety are constantly increasing, and structural analysis is becoming more refined. Therefore, researching more robust, accurate, and efficient stability analysis methods is urgently needed. Summary of the Invention

[0004] To address the aforementioned technical problems, this invention proposes a method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method, thereby resolving the issues existing in the prior art.

[0005] To achieve the above objectives, this invention provides a method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method, comprising:

[0006] S1. Initialize state variables and anchored secant predictor;

[0007] S2. Based on the state variables and the anchored secant predictor, predict the reference displacement, calibration displacement, and load increment factor;

[0008] S3. Calculate the displacement increment of the current iteration step based on the reference displacement, calibration displacement and load increment factor, and update the total structural displacement and total load.

[0009] S4. Update the structural geometry and element internal forces based on the displacement increment of the current iteration step;

[0010] S5. Based on the updated structural geometry and element internal forces, calculate the difference between the total internal forces of the structure and the external load to obtain the structural unbalanced force, and perform error verification based on the structural unbalanced force.

[0011] S6. Repeat S1-S5 until the unbalanced force of the structure meets the equilibrium condition, and start the next incremental step; when the incremental analysis is executed until the structural load or displacement meets the user's requirements, output the buckling path after multiple critical points.

[0012] Optionally, the state variables include: node coordinates, element coordinate basis, node displacement, element internal force, external load, load increment factor, and element stiffness matrix.

[0013] Optionally, the process of predicting the reference displacement, calibration displacement, and load increment factor includes:

[0014] Based on the state variables and the anchored secant predictor, the displacement increment of the current incremental step is predicted, and the initial displacement increment includes: a reference displacement and a calibration displacement.

[0015] The load increment factor for the current iteration step is calculated using the generalized displacement control method based on the reference displacement and the calibration displacement.

[0016] Optionally, the process of predicting the displacement increment for the current increment step includes:

[0017] If the current iteration step is 1, the reference displacement and calibration displacement are calculated based on the anchored secant predictor; wherein, the expression for calculating the reference displacement and calibration displacement based on the anchored secant predictor is:

[0018]

[0019] If the current iteration step is greater than 1, the reference displacement and calibration displacement are calculated based on the definition; wherein, the expressions for calculating the reference displacement and calibration displacement based on the definition are:

[0020]

[0021] In the formula, This represents the reference displacement when the current iteration step is 1. The calibrated displacement, {ΔU}, represents the calibration displacement when the current iteration step is 1. i-1} represents the displacement increment for each iteration in the (i-1)th increment step. The cumulative value, λ i-1 This represents the incremental load factor for each iteration in the (i-1)th incremental step. The cumulative value, This represents the reference displacement when the current iteration step is greater than 1. This represents the calibration displacement when the current iteration step is greater than 1. This represents the structural stiffness matrix at the beginning of the i-th increment step and the j-th iteration step. and These represent the user-preset reference load and the unbalanced force calculated in the previous iteration step, respectively.

[0022] Optionally, the expression for calculating the load increment factor of the current iteration step using the generalized displacement control method is as follows:

[0023] If the current iteration step is 1:

[0024]

[0025] In the formula, The initial load increment factor preset for the user; GSP is a generalized stiffness parameter.

[0026] If the current iteration step is greater than 1:

[0027]

[0028] In the formula, This represents the load increment factor. This represents the reference displacement of the initial increment step, the (i-1)th increment step, and the first iteration step.

[0029] Optionally, the process of calculating the displacement increment of the current iteration step and updating the total structural displacement and total load includes:

[0030] The displacement increment of the current iteration step is obtained by superimposing the displacement increment caused by the external load with the calibration displacement based on the superposition principle; wherein, the displacement increment caused by the external load is the product of the load increment factor and the reference displacement.

[0031] The external load increment is obtained by multiplying the load increment factor by the reference load.

[0032] The total structural displacement and total load are updated based on the external load increment and the displacement increment of the current iteration step.

[0033] Optionally, the expressions for updating the total structural displacement and total load based on the external load increment and the displacement increment of the current iteration step are as follows:

[0034]

[0035] In the formula, and Let be the total loads applied to the structure at the (j-1)th and jth iteration steps, respectively, for the i-th increment step. and These represent the total displacement responses of the structure from its initial state to the (j-1)th and jth iteration steps, respectively. Indicates the increment of external load. This represents the displacement increment of the current iteration step.

[0036] Optionally, the process of updating the structural geometry and element internal forces based on the displacement increment of the current iteration step includes:

[0037] The updated nodal coordinates and updated element coordinate basis are obtained by updating the structural geometry based on the updated total structural displacement and total load and the state variables;

[0038] Based on the updated node coordinates and the updated element coordinate base, the elastic stiffness matrix and geometric stiffness matrix of each element are calculated using the state variables to obtain the updated element elastic stiffness matrix and geometric stiffness matrix.

[0039] The updated element internal forces are obtained based on the updated node coordinates, the updated element coordinate base, the updated element elastic stiffness matrix, and the updated geometric stiffness matrix.

[0040] Optionally, the equilibrium condition is:

[0041]

[0042] In the formula, ε represents the unbalanced force in the structure. r Indicates the relative error limit. This represents the total load applied to the structure in the j-th iteration step.

[0043] Compared with the prior art, the present invention has the following advantages and technical effects:

[0044] The method of this invention is based on nonlinear finite element theory and establishes an efficient, stable and accurate method for structural nonlinearity and post-buckling analysis. It can effectively identify the post-buckling path and multiple critical points of the structure, which is of great significance for predicting the structural bearing capacity or ultimate strength. Attached Figure Description

[0045] The accompanying drawings, which form part of this application, are used to provide a further understanding of this application. The illustrative embodiments and descriptions of this application are used to explain this application and do not constitute an undue limitation of this application. In the drawings:

[0046] Figure 1 This is a flowchart of the method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method according to an embodiment of the present invention.

[0047] Figure 2 This embodiment of the invention is based on an incremental-iterative mechanism for updating the Lagrange formula;

[0048] Figure 3The prediction model proposed in this embodiment of the invention is: an anchored secant predictor;

[0049] Figure 4 This is a schematic diagram of a two-dimensional hinged deep circular arch structure in an embodiment of the present invention;

[0050] Figure 5 This refers to the post-buckling response under multiple critical points under symmetrical loading conditions calculated based on the method of this invention in this embodiment of the invention.

[0051] Figure 6 This refers to the post-buckling response under multiple critical points under eccentric loading conditions calculated based on the method of this invention in this embodiment of the invention.

[0052] Figure 7 This is a schematic diagram of a three-dimensional open spherical shell structure in an embodiment of the present invention;

[0053] Figure 8 This refers to the multiple critical point post-buckling response calculated based on the method of this invention in an embodiment of the invention. Detailed Implementation

[0054] It should be noted that, unless otherwise specified, the embodiments and features described in this application can be combined with each other. This application will now be described in detail with reference to the accompanying drawings and embodiments.

[0055] It should be noted that the steps shown in the flowchart in the accompanying drawings can be executed in a computer system such as a set of computer-executable instructions, and although a logical order is shown in the flowchart, in some cases the steps shown or described may be executed in a different order than that shown here.

[0056] Example 1

[0057] like Figure 1 As shown, this embodiment provides a method for rapidly and stably solving the post-buckling path of multiple critical points based on the generalized displacement control method. It is a solution method for structural nonlinearity and post-buckling problems based on the nonlinear finite element method, which is particularly beneficial for, but not limited to, using only elastic stiffness for stability analysis.

[0058] This invention novelly proposes an "anchored secant predictor" to predict the first iteration step of each increment step. That is, the secant stiffness of the previous increment step is used to approximate the tangent stiffness of the current first iteration step. This can complete the prediction stage of the first iteration step without the need for a stiffness matrix, which greatly improves the computational efficiency.

[0059] This invention combines the "anchored secant predictor" with the "generalized displacement control method" to achieve a fast, stable and accurate solution for buckling paths after multiple critical points.

[0060] This invention allows users to perform post-buckling analysis of structures using only elastic stiffness in moderately nonlinear problems, while overcoming the efficiency reduction problem that occurs when using only the elastic stiffness matrix for buckling analysis in the prior art.

[0061] In solving highly nonlinear problems, this invention strengthens the orthogonality condition of the generalized displacement control method by adding a qualified geometric stiffness matrix of the rigid body in subsequent iterations, so as to ensure that the iteration can converge stably even when the post-buckling path has a large curvature.

[0062] Aside from the necessary initialization definitions (including the initial settings of state variables and the initial conditions for the anchored secant predictor), the solution process can be mainly divided into three stages: prediction stage, calibration stage, and error verification stage. To maintain the generality of the method, the calculation process of the j-th iteration step in the i-th increment step is described here:

[0063] Prediction phase:

[0064] Step 1: Inherit the state variables from the previous iteration, including nodal coordinates, element coordinate basis, nodal displacements, element internal forces, external loads, initial load increment factors, and element stiffness matrix. If the current iteration j>1, further calculate the coordinate transformation matrix [T] based on the state variables, assemble the element stiffness matrix into the structural stiffness matrix [K], and simultaneously correct the structural stiffness matrix according to the boundary conditions.

[0065] Step 2: Predict the structural reference displacement and calibration displacement If the current iteration step j = 1, the reference displacement With calibration displacement The calculation will be performed using the proposed "anchored secant predictor":

[0066]

[0067] If the current iteration step j > 1, the reference displacement With calibration displacement Then, by definition, the following calculation is performed:

[0068]

[0069] In the above formula {ΔU i-1} represents the displacement increment for each iteration in the (i-1)th increment step. The cumulative value, λ i-1 The incremental load factor for each iteration in the (i-1)th incremental step The cumulative value. The structural stiffness matrix at the start of the i-th incremental step and the j-th iteration step is obtained by assembling the element stiffness matrix calculated in the (j-1)-th iteration step. and These represent the user-preset reference load and the unbalanced force calculated in the previous iteration step, respectively. (If the current iteration step j=1, the structure is in equilibrium, and there is an unbalanced force at this time.) Therefore, the calibration displacement calculated from the unbalanced force ).

[0070] Step 3: Using the calculated reference displacement and calibration displacement, calculate the load increment factor according to the generalized displacement control method. The load increment factor for the first iteration step is:

[0071]

[0072] In the above formula The initial load increment factor preset for the user determines the length of the increment step. GSP is a generalized stiffness parameter, defined as the ratio of the initial displacement increment to the current displacement increment of the structure.

[0073]

[0074] in and These are the reference displacements for the first iteration of the initial increment step, the first iteration of the (i-1)th increment step, and the first iteration of the ith increment step, respectively. For the j-th iteration step (j>1), the load increment factor is:

[0075]

[0076] and These are the reference displacement and calibration displacement for the j-th iteration step of the i-th increment step, respectively.

[0077] Step 4: Using the superposition principle, the displacement increment caused by the external load (expressed as the load increment factor) is calculated. With reference displacement The product of ( ) and the calibration displacement caused by the unbalanced force. By superimposing these values, the displacement increment of the i-th increment step and the j-th iteration step can be obtained. The corresponding external load increment Represented as load increment factor With reference load The product of, i.e. This updates the total structural displacement and total structural load in the current iteration step:

[0078]

[0079] In the above formula and These represent the total loads applied to the structure at the (j-1)th and jth iteration steps, respectively, in the i-th increment step; the corresponding... and These represent the total displacements generated by the structure from its initial state to the (j-1)th and (j)th iteration steps, respectively.

[0080] Calibration phase:

[0081] Step 5: Update the structural geometry (including nodal coordinates and element coordinate basis) based on the predicted structural displacement increments from Step 4. Based on the updated geometry, calculate the coordinate transformation matrix [T] for each element sequentially, and utilize the displacement increments obtained in Step 4. The displacement increment of each element is obtained through the coordinate transformation matrix [T].

[0082] Step 6: Traverse each element and calculate the elastic stiffness matrix [k] of the element under each updated geometry. e ] and geometric stiffness matrix [k g (For moderately nonlinear problems, the geometric stiffness matrix can be omitted. However, when solving highly nonlinear problems with large curvature at the springback point in the post-buckling path, a qualified geometric stiffness matrix for the rigid body must be included. Here, for the sake of program versatility, the geometric stiffness matrix [k] will still be calculated.) g (The steps are listed below)

[0083] Step 7: Based on the rigid body criterion and the small deformation increment theory, traverse all elements and update the element internal forces under the deformed configuration:

[0084]

[0085] In the above formula The internal forces of the elements inherited from step 1 in the (j-1)th iteration before the element geometry was updated, and are referenced to the element coordinate base of the element with the updated geometry; while The element internal forces after updating the structural geometry in step 5 are referenced to the element coordinate base after the structural geometry update. (Here, symbols...) The superscript "1" in the middle left indicates that the variable belongs to the (h-1)th iteration step, while the subscript "1" in the left left indicates that the variable is referenced to the coordinate base of the (j-1)th iteration step; and the symbol... Similarly, (representing the j-th iteration step)

[0086] Error verification stage:

[0087] Step 8: Traverse all units and update the internal forces of each unit. Stored for later use when updating the element internal forces in the next iteration step, i.e., as the element internal forces in the next iteration step. Subsequently, the internal forces of each unit will be... The total internal forces of the structure are obtained by transforming the coordinates to the structural coordinate base using the coordinate transformation matrix [T] and then superimposing the internal forces of all elements.

[0088] Step 9: Calculate the unbalanced forces of the structure based on the difference between the total external load and the total internal force. And determine whether the structure has reached equilibrium. (typical relative error limit ε) r (Can be set to 0.0001) If the balance condition is not met, return to step 1 to start the (j+1)th iteration step; if the balance condition is met, the i-th increment step is considered complete, and return to step 1 to start the (i+1)-th increment step 1 to start the first iteration step.

[0089] Furthermore, the theoretical basis of the technical solution includes:

[0090] In continuity mechanics, the motion and deformation of a structure can be described by three configurations, with the initial configuration C being the first configuration. 0 The configuration C is known from the previous step. i-1 and the current unknown configuration C i Based on the updated Lagrange formula, select the known configuration C from the previous step. i-1 Using the reference configuration to describe the structure from configuration C i-1 Position C i Incremental behavior between them. Figure 3 This demonstrates an incremental-iterative mechanism based on the Lagrange formula, such as... Figure 2 As shown, the theoretical basis of this invention (in the figure and subsequent theory, the superscript of each variable indicates the increment step when it appears, and the subscript indicates the iteration step):

[0091] (1) The structure is in the known configuration C i-1 At time Figure 3 Point a), under load {P} i-1 Under the action of}, it is in equilibrium, and its corresponding displacement response is {U i-1 Apply a small load increment {ΔP} to the structure. i} Make the load from {P i-1} changes to {P i}

[0092] (2) Prediction stage: based on stiffness equation The structure can be solved under load increment {ΔP} i The displacement increment {ΔU} under} i (At this point, the trajectory ab is reached) Figure 3 (Point b), and simultaneously calculate the displacement increment {Δu} of each element in the element coordinate system through coordinate transformation. i}

[0093] (3) Calibration stage: Based on the required displacement increment {ΔU} i Update structural geometry (from) Figure 3 From point a to point d), and calculate the total internal force of the element (from point a to point d). Figure 3 (Represented by line segment ce). The calculation of the total internal force of the element includes two parts: the change in internal force of the element caused by the rotation of the rigid body (by... Figure 3 (from line segment oa to line segment cd) and the element force increment caused by natural deformation (by Figure 3 (represented by the line segment de).

[0094] (4) Error verification stage: The total internal force {F} of the structure is obtained by calculating and superimposing the internal forces of each element based on coordinate transformation. i Then, by calculating the total external load and the total internal forces of the elements, the unbalanced forces of the structure are obtained, i.e., {R}. i}={P i}-{F i}(Depend on Figure 3 (represented by line segment eb)

[0095] (5) Next iteration: Apply the unbalanced force as an external load to the structure and start a new iteration step to calculate the structural displacement and element force; if the unbalanced force meets the preset allowable conditions, the structure is considered to be in equilibrium, the calculation of the i-th increment step is terminated, and the next i+1-th increment step begins.

[0096] Based on the prediction stage in the above incremental-iterative mechanism, the incremental stiffness equation for the i-th incremental step and the j-th iteration step can be expressed in general form:

[0097]

[0098] in It is the structural stiffness matrix starting from the i-th increment step and the j-th iteration step. It is the displacement increment to be determined. It is the set reference load. Let be the load increment factor to be determined. This is the unbalanced force of the structure in the (j-1)th iteration step. Typically, the displacement increment to be determined can be divided into... and Two parts. This is called the reference displacement, determined by the reference load. Caused by; and This is called the calibration displacement, caused by unbalanced forces. Caused by:

[0099]

[0100] Equation (9) contains a total of N equations, but it contains N+1 unknowns (including the displacement increments of N degrees of freedom). With constant Therefore, in order to solve equation (9), a constraint equation needs to be added, the generalized form of which is:

[0101]

[0102] In the above formula, k, {C} and These are all constraint parameters that need to be assumed in advance; different constraint parameters can distinguish different path tracing methods. In the generalized displacement control method, these are set... and The initial conditions are subsequently determined. Combining equations (9) and (13), the load increment factor can be obtained as follows:

[0103]

[0104] At the start of the first iteration of each incremental step, the structure is in equilibrium, and the displacement is calibrated at this time. Equation (14) is transformed into:

[0105]

[0106] In the generalized displacement control method, it is assumed that... It remains constant in each increment step, that is, in the first iteration step j=1. In subsequent iterations j≥2, Therefore, for the second iteration step:

[0107]

[0108] Substitute initial conditions Therefore, the solution is obtained. Substituting it back into equation (15) yields:

[0109]

[0110] GSP stands for Generalized Stiffness Parameter, defined as the ratio of the initial displacement increment to the current displacement increment. It is used in the generalized displacement control method to reflect changes in structural stiffness.

[0111]

[0112] The generalized stiffness parameter (GSP), used as the control parameter in the generalized displacement control method, is characterized by being negative only in the first increment step after passing the limit point. At this point, the direction of the load increment should be opposite to that of the previous increment step. Simply put, when GSP > 0 (at which point the reference displacement of the (i-1)th increment step is...), With the reference displacement of step i (At an acute angle), the load increment factor for the i-th increment step The sign of the load increment should be the same as that of the (i-1)th step, except when GSP < 0 is detected (in which case the reference displacement of the (i-1)th increment step is...). With the reference displacement of step i (The angle is obtuse), so the direction of the load increment factor should be opposite to the sign of the load increment factor in step i-1.

[0113] Traditionally, for the prediction phase of each iteration step, the reference displacement is... and calibration displacement These can be obtained from the definitions (10) and (11), respectively. Substituting them into equations (17) and (16) yields the load increment factors, respectively. (Iteration step 1) and (In the j-th iteration step, j>1). However, for the 1st iteration step, if the reference displacement... The proposed anchored secant predictor can significantly reduce computation time while ensuring prediction accuracy.

[0114] Furthermore, the principles of anchored secant prediction include:

[0115] Extending the concept of secant stiffness from the previous increment step to predict the first iteration step of the i-th increment step, as follows: Figure 3 As shown. For loads from {P} i-1} to {P i In the i-th increment step, line segments ac and ab are the first iteration steps calculated using secant stiffness and tangent stiffness, respectively. Their trajectories will pass through subsequent iteration points m and f, and eventually converge to the vicinity of point s. When the increment step size is relatively small, i.e., point n is very close to point a, the secant stiffness anchored in the (i-1)-th increment step (the slope of line segment na or ac) will also be very close to the tangent stiffness of the first iteration step in the current i-th increment step (the slope of line segment ab). It is worth noting that, Figure 5 Points c and b in the diagram are very close, but for ease of explanation, they are exaggerated here. Based on Figure 5 The secant stiffness equation for Δnar at the (i-1)th increment step can be written as:

[0116]

[0117] In the above formula It is the secant stiffness at the (i-1)th increment step; {ΔU i-1} and λ i-1 These are the cumulative displacement increment and cumulative load increment factor for the (i-1)th increment step, respectively. By definition, the cumulative load {ΔP}i-1} can be represented as a reference load. With the cumulative load increment factor λ i-1 The product of.

[0118] Based on the reference displacement of the first iteration step The calculation is defined by equation (10). Taking advantage of the fact that the secant stiffness and tangent stiffness are approximately equal in the sense of incremental step, the calculation of the reference displacement in the first iteration step can be written as:

[0119]

[0120] Combining equations (19) and (20), we can simplify to obtain:

[0121]

[0122] Equation (21) is called the anchored secant predictor. It is worth noting that in order to start the anchored secant predictor, the first iteration of the first increment step should be calculated using the definition (10). The characteristic of the anchored secant predictor is that once the (i-1)th increment step is completed, the first iteration of the i-th increment step can be predicted directly using the concept of secant stiffness. Note that Equation (21) does not require any information related to the stiffness equation; it only needs to perform assignment operations using the data from the previous increment step. Therefore, the steps of assembling the structural stiffness matrix and solving the linear equation system in the first iteration step can be omitted, thus greatly saving computation time.

[0123] Furthermore, the orthogonal properties of the generalized displacement control method include:

[0124] The purpose of subsequent iterations is to eliminate the errors caused by the approximation of the incremental method in the first iteration. Based on the generalized displacement control method, the constraint equation for the j-th iteration (j>1) is:

[0125]

[0126] In the above formula It is the undetermined displacement increment in the current j-th iteration step. This is the reference displacement increment of the first iteration step in the (i-1)th increment step, which can be regarded as the guiding vector of the iteration process. The constraint equation (22) can be regarded as the displacement to be located. With guidance vector The orthogonality condition between them, i.e., subsequent iterations and Always with the guiding displacement The directions are perpendicular.

[0127] Furthermore, element internal force calibration based on rigid body criteria includes:

[0128] Element internal force calibration at the element level mainly involves two aspects: (1) recovering the element force increment {Δf} from the displacement increment; (2) updating the total element internal force acting on the element. After the internal force calibration phase is completed, the total internal force of the element is used to calculate the structural unbalanced force {R}. It is conceivable that if the total internal force of the element... Inaccurate updates will naturally result in ill-conditioned calculated unbalanced forces; subsequent iterations will therefore fail to converge to the correct position. Thus, the update of the total internal forces of the element during the element internal force calibration stage determines the accuracy of the obtained solution.

[0129] Based on physical objectivity, the displacement {Δu} of any finite element can be divided into two parts: rigid body displacement {Δu}. r and natural deformation {Δu} n When the structure becomes unstable, the rigid body of the element rotates {Δu}. r The rigid body displacement accounts for the vast majority of the element displacement increment {Δu}; in contrast, the natural deformation of the element is much smaller. This objective phenomenon requires that, in incremental-iterative analysis, the change in element internal force caused by rigid body displacement must be accurately considered, as it represents the main part of the element displacement increment; conversely, the element force increment caused by natural deformation can be regarded as the element force increment within a small deformation range, and can be calculated using only elastic stiffness.

[0130] {Δf}=[k e ]{Δu} n (twenty three)

[0131] Since the elastic stiffness matrix needs to meet the requirements of piecewise testing, no element force increment will be generated under rigid body displacement. Therefore, the above equation (23) can also be rewritten as:

[0132] {Δf}=[k e ]{Δu} n =[k e ]({Δu}-{Δu} r )=[k e ]{Δu} (24)

[0133] On the other hand, the changes in internal forces of an element caused by rigid body displacement can be accurately handled according to the rigid body criterion. The rigid body criterion is a simple and intuitive physical principle, universally effective for testing the mass of a finite element subjected to initial forces. For a finite element subjected to a set of balanced forces before rigid body rotation, the internal forces acting on the element during the rotation will rotate along with the rigid body. After the rotation, the element remains in equilibrium, and the forces acting on the element have changed only in direction, not magnitude, compared to before the rotation. Mathematically, the rigid body criterion can be expressed as:

[0134]

[0135] The changes in element internal forces caused by rigid body rotation, as indicated by equation (25), can be completely and accurately considered. The procedure is also very simple, requiring no additional coordinate transformations; only the element internal forces from the previous iteration step need to be considered. Numerically, it can be directly regarded as the internal force of the element after the rigid body rotates. That's it. By utilizing the changes in element internal forces caused by superimposed rigid body displacements and the increments in element internal forces caused by natural deformation, the total internal forces of the elements in the updated geometry can be obtained:

[0136]

[0137] The advantage of the above formula is that it only requires the elastic stiffness matrix to complete the internal force update of the element, without being affected by the approximation of the geometric stiffness derivation.

[0138] The method of this invention is based on nonlinear finite element theory and establishes an efficient, stable and accurate method for structural nonlinearity and post-buckling analysis. It can effectively identify the post-buckling path and multiple critical points of the structure, which is of great significance for predicting the structural bearing capacity or ultimate strength.

[0139] Example 2

[0140] As a specific implementation of this embodiment, further numerical verification of the technical solution of the present invention is performed.

[0141] Example 1 parameters:

[0142] The geometry of a two-dimensional hinged deep circular arch structure is as follows: Figure 4 As shown, under the action of a vertical load P, the arch span L = 254 cm, and the interface moment of inertia Iz = 41.62 cm. 4 The cross-sectional area A = 64.52 cm² 2 The elastic modulus E = 1378 kPa. In the finite element analysis, the deep circular arch is first divided into 25 equal two-dimensional beam elements; then, the arch crown element is divided into two elements to provide a central vertex. Therefore, a total of 26 two-dimensional beam elements are used to simulate this problem. In this embodiment, two loading conditions will be considered: one is symmetrical loading, where the vertical load P will be applied at the central node of the arch crown; and the other is eccentric loading, where the vertical load P will be applied to adjacent nodes that are off-center from the central node. A reference load is set in the incremental analysis. Initial load increment factor It should be added that the numerical simulations of all embodiments in this invention are performed in 13th Gen. Core TMThe i5-13600KF is equipped with a 3.50GHz CPU and 32.0GB of RAM, running on Windows 11, 64-bit operating system.

[0143] To verify the improved computational efficiency of the anchored secant predictor proposed in this invention, we performed numerical simulations and compared the computational efficiency with that of the traditional predictor. The traditional predictor is defined in equation (10). Figure 5 (a) shows the load-deflection curves of two predicted deep circular arches under symmetrical load conditions. Figure 5 (b) shows the GSP-deflection curves for the two predictors under symmetrical loading conditions. Figure 6 (a) shows the load-deflection curves of two predicted deep circular arches under eccentric loading conditions. Figure 6 (b) shows the GSP-deflection curves of the two predictors under symmetrical load conditions.

[0144] from Figure 5 (a) and Figure 6 (a) It can be seen that both predictors can effectively track the post-buckling path of the structure and identify multiple critical points; in addition... Figure 5 (b) and Figure 6 (b) also shows that the two predictors perform essentially the same in predicting the generalized stiffness parameter GSP, and both can identify the zeros on the GSP curve (which correspond to the extreme points in the post-buckling path). The computational efficiency comparison between the anchored secant predictor and the traditional predictor under the two loading conditions is shown in Tables 1 and 2:

[0145] Table 1

[0146]

[0147]

[0148] Table 2

[0149]

[0150] The calculation results show that, under both central and eccentric load conditions, the incremental steps and average iteration steps required for the two predictors are comparable. However, it is evident that the calculation time of the anchored secant predictor is shorter than that of the traditional predictor, demonstrating the advantage of the proposed anchored secant predictor in improving computational efficiency.

[0151] Example 2 parameters:

[0152] The geometry of a three-dimensional open spherical shell structure is as follows: Figure 7As shown, the shell edges are all hinged, and a vertical load P acts at the center. The shell has a side length of 2a = 1569.8 mm, a thickness of t = 99.45 mm, a radius of curvature R = 2540 mm, a Poisson's ratio v = 0.3, and an elastic modulus E = 6895 kPa. Due to the axisymmetric nature of the shell, only one-quarter of the shell is modeled as an 8×8 mesh for analysis in the finite element modeling. A reference load is set in the incremental analysis. Initial load increment factor

[0153] To verify the improvement in computational efficiency of the method of this invention, we performed numerical simulations and compared the computational efficiency with that of two traditional predictors. Both traditional predictors are defined by equation (10), where the stiffness matrix of traditional predictor 1 is... It also includes the elastic stiffness matrix [K] e ] and geometric stiffness matrix [K g ]; while the traditional predictor 2 only contains the elastic stiffness matrix [K e ]. Figure 8 (a) shows the load-deflection curves of the three predicted sub-spherical shells under the central load condition. Figure 8 (b) shows the GSP-deflection curves of the three predicted subspherical shells under the central load condition.

[0154] from Figure 8 (a) It can be seen that all three predictors can effectively track the post-buckling path of the structure and identify multiple critical points; however, it should be noted that Figure 8 (b) The difference between the GSP curve of traditional predictor 2 and the other two is not negligible. This is because a large error occurred when predicting the displacement of the first iteration step using only the predictor with elastic stiffness, leading to errors in the GSP calculation. This also explains why using only elastic stiffness for prediction results in reduced computational efficiency, as the error in the first iteration step requires more subsequent iterations to eliminate, thus wasting computation. However, it can be seen that the anchored secant predictor proposed in this invention can cleverly solve the efficiency reduction problem that occurs when using only elastic stiffness. Furthermore, the computational efficiency results of the three predictors are shown in Table 3:

[0155] Table 3

[0156]

[0157] The efficiency improvement of the anchored secant predictor in Example 2 is even more significant, with computational efficiency nearly twice that of the conventional predictor. Furthermore, by utilizing the anchored secant predictor, the problem of low computational efficiency reported in the literature can be overcome when using only elastic stiffness.

[0158] The above are merely preferred embodiments of this application, but the scope of protection of this application is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in this application should be included within the scope of protection of this application. Therefore, the scope of protection of this application should be determined by the scope of the claims.

Claims

1. A method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method, characterized in that, Includes the following steps: S1. Initialize state variables and anchored secant predictor; S2. Based on the state variables and the anchored secant predictor, predict the reference displacement, calibration displacement, and load increment factor; S3. Calculate the displacement increment of the current iteration step based on the reference displacement, calibration displacement and load increment factor, and update the total structural displacement and total load. S4. Update the structural geometry and element internal forces based on the displacement increment of the current iteration step; S5. Based on the updated structural geometry and element internal forces, calculate the difference between the total internal forces of the structure and the external load to obtain the structural unbalanced force, and perform error verification based on the structural unbalanced force. S6. Repeat S1-S5 until the unbalanced force of the structure meets the equilibrium condition, and start the next incremental step. The incremental analysis is executed until the structural load or displacement meets the user's requirements, and then the buckling path after multiple critical points is output. The process of predicting reference displacement, calibrating displacement, and load increment factor includes: Based on the state variables and the anchored secant predictor, the displacement increment of the current increment step is predicted. The initial displacement increment includes: reference displacement and calibration displacement. The load increment factor for the current iteration step is calculated using the generalized displacement control method based on the reference displacement and the calibration displacement. The process of predicting the displacement increment of the current increment step includes: If the current iteration step is 1, the reference displacement and calibration displacement are calculated based on the anchored secant predictor; wherein, the expression for calculating the reference displacement and calibration displacement based on the anchored secant predictor is: If the current iteration step is greater than 1, the reference displacement and calibration displacement are calculated based on the definition; wherein, the expressions for calculating the reference displacement and calibration displacement based on the definition are: In the formula, This represents the reference displacement when the current iteration step is 1. This represents the calibration displacement when the current iteration step is 1. Indicates the first Incremental displacement in each iteration of the incremental step The cumulative value, Indicates the first Incremental load factor for each iteration of the incremental step The cumulative value, This represents the reference displacement when the current iteration step is greater than 1. This represents the calibration displacement when the current iteration step is greater than 1. Indicates the first i Incremental step number j The structural stiffness matrix at the start of the iteration step. and These represent the user-preset reference load and the unbalanced force calculated in the previous iteration step, respectively. The expression for calculating the load increment factor of the current iteration step using the generalized displacement control method is as follows: If the current iteration step is 1: In the formula, The initial load increment factor preset for the user; GSP is a generalized stiffness parameter. If the current iteration step is greater than 1: In the formula, Indicates the load increment factor. Indicates the initial increment step. The reference displacement for the first iteration of the increment step; The process of calculating the displacement increment at the current iteration step and updating the total structural displacement and total load includes: The displacement increment of the current iteration step is obtained by superimposing the displacement increment caused by the external load with the calibration displacement based on the superposition principle; wherein, the displacement increment caused by the external load is the product of the load increment factor and the reference displacement. The external load increment is obtained by multiplying the load increment factor by the reference load. The total structural displacement and total load are updated based on the external load increment and the displacement increment of the current iteration step.

2. The method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method according to claim 1, characterized in that, The state variables include: node coordinates, element coordinate basis, node displacement, element internal force, external load, initial load increment factor, and element stiffness matrix.

3. The method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method according to claim 1, characterized in that, The expressions for updating the total structural displacement and total load based on the external load increment and the displacement increment of the current iteration step are as follows: In the formula, and The first The first incremental step Iteration step and the The total load applied to the structure in each iteration step and The structure is from its initial state to the... Incremental step number Iteration step and the The total displacement response generated by the iterative step structure. Indicates the increment of external load. This represents the displacement increment of the current iteration step.

4. The method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method according to claim 1, characterized in that, The process of updating the structural geometry and element internal forces based on the displacement increment of the current iteration step includes: The updated nodal coordinates and updated element coordinate basis are obtained by updating the structural geometry based on the updated total structural displacement and total load and the state variables; Based on the updated node coordinates and the updated element coordinate base, the elastic stiffness matrix and geometric stiffness matrix of each element are calculated using the state variables to obtain the updated element elastic stiffness matrix and geometric stiffness matrix. The updated element internal forces are obtained based on the updated node coordinates, the updated element coordinate base, the updated element elastic stiffness matrix, and the updated geometric stiffness matrix.

5. The method for rapidly and stably solving the buckling path after multiple critical points based on the generalized displacement control method according to claim 1, characterized in that, The equilibrium condition is: In the formula, Indicates unbalanced forces in the structure. Indicates the relative error limit. Indicates the first The total load applied to the structure in each iteration step.

Images