Topological Optimization Method and System Considering Stability under Design-Related Loads

Abstract

This invention provides a topology optimization method and system that considers stability under design-related loads. This relates to the technical field of structural topology optimization. By using shape derivative theory and augmented Lagrangian operators, a velocity field is obtained that takes structural stability into account under design-related loads. Furthermore, a regularization term is introduced into the H-J equation to control excessive structural differentiation during the optimization process. This method avoids excessive differentiation of fine structural branches within the structure during the optimization process, thereby improving the rationality and manufacturability of the final structure.

Description

Technical Field

[0001] The present invention relates to the technical field of structural topology optimization, and in particular to a topology optimization method that considers structural stability under the action of design-related loads, and in particular to a topology optimization method and system that considers stability under the action of design-related loads. Background Art

[0002] In the optimization design of deep-sea pressure-resistant structures, the loads applied to the structure are different from the fixed loads used in traditional topology optimization design. This type of hydrostatic pressure changes as the structural contour changes. Therefore, this type of load is also called a design-related load. In deep-sea environments, pressure-resistant structures are subject to hydrostatic pressure, and the failure forms of the structure are generally the following two: material failure and structural instability. Structural instability is particularly critical for the design of deep-sea pressure-resistant structures. Structural instability, also known as structural buckling, refers to the phenomenon that when the external force increases to a certain amount, the structure begins to lose its stable equilibrium state. With a slight disturbance, the structural deformation will increase rapidly, causing the structure to fail. Structural instability can be divided into two categories: bifurcation point instability and extreme point instability. The difference between the two is that the former is a mathematical instability that does not consider the mechanical properties of the structure, while the latter takes into account the changes in material properties during the stress process. While the analysis results obtained using extreme point instability analysis are more accurate than those obtained using bifurcation point instability theory, the results obtained using bifurcation point instability theory are computationally simple and highly efficient. Therefore, this method is widely used in the field of topology optimization design. Topology optimization is a mathematical method that optimizes the distribution of materials within a given area based on given load conditions, constraints, and performance indicators. Using the shape of the structure as a design variable, this optimization method seeks to achieve an optimal topological form that meets design requirements. This method significantly changes the mindset of engineering designers and provides guidance for engineering practice.

[0003] Unlike traditional topology optimization problems involving fixed loads, design-related loads vary with the structural contour. These problems involve searching for the structural boundary contours to find the optimal external configuration. Therefore, these problems play an important role in the engineering application of deep-sea pressure-resistant structures.

[0004] In the field of topology optimization, the goal of topology optimization under design-related loads is to find the optimal external and internal structure for structures subjected to hydrostatic pressure. Existing topology optimization methods that consider stability under design-related loads are based on density methods. However, because the structural boundaries obtained by density methods are too vague, a significant amount of time is spent on boundary search during the computational process, reducing computational efficiency.

[0005] Finding a topology optimization method that can provide a clear structural boundary at all times during the optimization process, thereby improving the computational efficiency of the stability considerations under design-related loads, is an urgent problem that needs to be solved. Summary of the Invention

[0006] In view of the defects in the prior art, the present invention provides a topology optimization method and system that considers stability under the action of design-related loads.

[0007] According to the present invention, a topology optimization method and system considering stability under design-related loads are provided, and the solution is as follows:

[0008] In a first aspect, a topology optimization method considering stability under design-related loads is provided, the method comprising:

[0009] Step S1: Initialize the level set function, input a set of radial basis functions and interpolate the level set function to obtain the interpolation coefficients of this set of radial basis functions;

[0010] Step S2: construct a topology optimization model with the structural shape as the design variable, the minimization of the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure as the optimization goal, and the volume of the material as the optimization constraint;

[0011] Step S3: According to the current level set function, search for the current force boundary contour of the structure, and use it to calculate the equivalent nodal force currently acting on the structure;

[0012] Step S4: Calculate the current objective function, shape derivatives of optimization constraints, and Lagrangian operators, solve the velocity field based on the calculation results to evolve and obtain a new level set function, and update the interpolation coefficients before the radial basis function;

[0013] Step S5: Determine the convergence. If the convergence condition is met, output the optimization results under the design-related loads with consideration of flexibility minimization and proceed to the next stage of optimization. Otherwise, return to step S3.

[0014] Step S6: Using the optimization result obtained in step S5 as the initial design, a topology optimization problem is performed with the goal of minimizing the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor and the material volume fraction as the optimization constraint. Based on the current level set function, the current force boundary contour of the structure is searched, and the equivalent nodal force currently acting on the structure is calculated.

[0015] Step S7: Under the current level set function and the load, calculate the inverse of the first four-order eigenvalue buckling factors of the structure;

[0016] Step S8: Calculate the current objective function, shape derivatives of optimization constraints, and Lagrangian operators, solve the velocity field based on the calculation results to evolve and obtain a new level set function, and update the interpolation coefficients before the radial basis function;

[0017] Step S9: Determine the convergence. If the convergence condition is met, output the stability maximization topology optimization result under the design-related load. Otherwise, return to step S6.

[0018] Preferably, in step S2, the expression of the topology optimization model for minimizing the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure includes:

[0019]

[0020]

[0021]

[0022]

[0023]

[0024]

[0025] Where J(Ω) is the current objective function value; u is the displacement field of the structure, p is the pressure on the Neumann boundary, p0 is the pressure value; v is the test function field; U ad is the displacement field space; J KS (μ i (Ω)) is the current objective function value; ρ is the relative volume density under the imaginary density model; n is the normal vector on the structure contour, Φ is the level set function, V(Ω)=∫ Ω dx represents the volume of the structure, and μ represents the inverse of the eigenvalue buckling factor of the structure. i (Ω) represents the inverse of the i-th order eigenvalue buckling factor of the structure.

[0026] Preferably, the characteristic buckling factor λ of the structure is obtained by the following formula:

[0027] ∫ Ω Ae(ψ)e(ψ)dx+λ∫ Ω Ae(u):q(ψ,ψ)dx=0

[0028] Where q(ψ,ψ)=ψ k,i ψ k,j , ψ is the mode of the structure after eigenvalue buckling occurs. Specifically, ψ k,i ψ k,jrepresents the inner product sum operation between the derivatives of each side of ψ; in order to simplify the optimization problem into a minimization problem, the eigenvalue buckling factor λ is taken as its reciprocal μ, and the reciprocal eigenvalue buckling factor μ of the structure is obtained by the following formula:

[0029] ∫ Ω Ae(u):q(ψ,ψ)dx+μ∫ Ω Ae(ψ)e(ψ)dx=0.

[0030] Preferably, in step S4, in order to smooth the optimization process and the obtained results and ensure the smooth progress of the optimization process, the inverse of the first four-order buckling factors is subjected to KS aggregation and used as the objective function:

[0031]

[0032] Among them, ρ is a parameter in the KS aggregation equation, which is used to make the result of KS aggregation approach the inverse of the first-order buckling factor μ1.

[0033] Preferably, in step S4, the evolution of the level set function is driven based on the improved HJ equation, and the specific formula is as follows:

[0034]

[0035] in, here Refers to a set of known CSRBF radial basis functions, α(t) represents the The interpolation coefficients, is a regularization term used to control the excessive differentiation of the internal structure during the optimization process, v n is the velocity field, which controls the deformation of the structure during the optimization process.

[0036] Preferably, the velocity field v is solved based on the shape derivative theory n , for the topology optimization problem of minimizing the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure, the Lagrangian equations of each model are listed to solve the adjoint fields of each problem:

[0037]

[0038]

[0039] Among them, the topology optimization problem of minimizing structural flexibility is a self-adjoint problem. For the topology optimization problem of minimizing the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factors of the structure, the adjoint variable χ of the inverse of the eigenvalue buckling factor is obtained by the following formula:

[0040]

[0041] in, Used to normalize the eigenvalue vector, is the displacement field under the current problem; e(χ) is the strain tensor accompanying the field χ.

[0042] Preferably, the shape derivative of the Lagrange equation of each model is:

[0043]

[0044]

[0045] Where div(*) represents the divergence of *; θ represents the descending direction of the shape derivative.

[0046] Preferably, the material volume fraction of the structure is controlled to reach the optimal constraint level by using the augmented Lagrangian operator, and the augmented Lagrangian operator is obtained by the following formula:

[0047]

[0048] Where V1 is the material volume fraction of the initial structure, μ is a fixed parameter, and V0 is the volume fraction constraint value; n r is the number of iterations of the Lagrangian operator form transformation, Λ k is the Lagrangian operator value at the kth iteration step, γ k Represents the value of γ in the γth iteration step. The role of γ is to control the material volume fraction of the structure to approach the volume fraction constraint value. k The optimization process continues to change as shown by the following formula:

[0049] γ k+1 =min(γ k +Δγ,γ max ),k>n r

[0050] Among them, Δγ is the increment of γ in the optimization process, γ max is the upper limit of γ.

[0051] Preferably, the velocity field v of each problem n Obtained by the following formula:

[0052]

[0053]

[0054] in,

[0055]

[0056] Where δ(φ) is the Dirac function obtained by the following formula:

[0057]

[0058] Among them, Δ is greater than 5, and Ω1 and Ω2 represent the loaded and unloaded areas respectively;

[0059] An implicit boundary search method is used to determine the load application boundary of the structural design. The search strategy is: perform directional scanning row by row or column by column along the direction of the external load, and perform lateral scanning on each node in the direction perpendicular to the scanning direction in turn until the node φ ≥ 0, and stop scanning. The node is marked as the main mark point, and then the adjacent unit nodes are marked as secondary mark points. Only at the main and secondary mark points, the hydrostatic pressure coefficient p0(x) is not 0; it is used to distinguish the areas with and without load Ω1 and Ω2. The equivalent nodal force of the loaded area is finally determined by the following formula:

[0060]

[0061] Among them, N i (ξ,η) is the shape function, ξ and η are the horizontal and vertical coordinates in the isoparametric coordinate system, and their values range from -1 to 1.

[0062] In a second aspect, a topology optimization system considering stability under design-related loads is provided, the system comprising:

[0063] Module M1: Initialize the level set function, input a set of radial basis functions and interpolate the level set function to obtain the interpolation coefficients of this set of radial basis functions;

[0064] Module M2: Construct a topology optimization model with the structural shape as the design variable, the minimization of the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure as the optimization objective, and the volume of the material as the optimization constraint;

[0065] Module M3: Based on the current level set function, search for the current force boundary contour of the structure and use it to calculate the equivalent nodal force currently acting on the structure;

[0066] Module M4: Calculates the current objective function, shape derivatives of optimization constraints, and Lagrangian operators. Based on the calculation results, it solves the velocity field to evolve and obtain a new level set function, and updates the interpolation coefficients before the radial basis function.

[0067] Module M5: Determine the convergence. If the convergence condition is met, the optimization results under the design-related loads with consideration of flexibility minimization are output and the next stage of optimization is carried out. Otherwise, return to module M3.

[0068] Module M6: Using the optimization results obtained in Module M5 as the initial design, the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor is minimized as the goal, and the material volume fraction is used as the optimization constraint for the topology optimization problem. Based on the current level set function, the current force boundary contour of the structure is searched, and the equivalent nodal force currently acting on the structure is calculated.

[0069] Module M7: Calculate the inverse of the first four eigenvalue buckling factors of the structure under the current level set function and load.

[0070] Module M8: Calculates the current objective function, shape derivatives of optimization constraints, and Lagrangian operators. Based on the calculation results, it solves the velocity field to evolve and obtain a new level set function, and updates the interpolation coefficients before the radial basis function.

[0071] Module M9: Determine the convergence. If the convergence conditions are met, output the stability maximization topology optimization results under the design-related loads. Otherwise, return to module M6.

[0072] On the third aspect, a method for preparing a pressure-resistant structure is provided, wherein a pressure-resistant structure corresponding to a topology optimization result of maximizing stability under the action of design-related loads output by the topology optimization method considering stability under the action of design-related loads is prepared.

[0073] Fourthly, a pressure-resistant structure is provided, which is prepared by outputting a pressure-resistant structure corresponding to a topology optimization result that maximizes stability under the action of design-related loads according to the topology optimization method that considers stability under the action of design-related loads.

[0074] Compared with the prior art, the present invention has the following beneficial effects:

[0075] 1. This invention addresses the problem of topology optimization under design-related loads by combining a parameterized level set method with stability issues to produce an optimized structure with smooth boundaries. Furthermore, since this method addresses the problem of topology optimization under design-related loads, it conforms to the model of hydrostatic pressure in the optimized design of deep-sea pressure-resistant structures and can meet the structural performance design requirements in deep-sea pressure-resistant structure manufacturing.

[0076] 2. This invention introduces a regularization term to control excessive structural differentiation during the optimization process. By controlling the perimeter of the structural contour, the regularization term avoids excessive differentiation of fine structural branches within the structure during the optimization process, thereby improving the rationality and manufacturability of the final structure.

[0077] Other beneficial effects of the present invention will be explained through the introduction of specific technical features and technical solutions in the specific implementation methods. Those skilled in the art should be able to understand the beneficial technical effects brought about by the introduction of these technical features and technical solutions. BRIEF DESCRIPTION OF THE DRAWINGS

[0078] Other features, objects and advantages of the present invention will become more apparent upon reading the detailed description of non-limiting embodiments with reference to the following drawings:

[0079] Figure 1 A schematic flow chart of a topology optimization method considering structural stability under design-related loads provided by the present invention;

[0080] Figure 2 This is a classification diagram of element types for the hypothetical material model in the level set method;

[0081] Figures 3(a) and 3(b) are schematic diagrams showing implementation of the load application boundary search for structural design;

[0082] Figure 4 The design area and boundary conditions of the water-pressure arch structure of the present invention;

[0083] FIG5(a) and FIG5(b) are diagrams showing the optimization results of the water-pressure-treated arched structure according to an embodiment of the present invention;

[0084] Figures 6(a) and 6(b) are first-order buckling mode diagrams of the optimized structure of the example of the present invention;

[0085] Figure 7 It is the iterative variation diagram of the first-order buckling factor of the structure;

[0086] Figure 8 It is the iterative change diagram of the structural flexibility;

[0087] Figure 9 Iterative variation diagram of volume fraction of structural materials. DETAILED DESCRIPTION

[0088] The present invention will be described in detail below with reference to specific embodiments. The following examples will help those skilled in the art to further understand the present invention, but are not intended to limit the present invention in any form. It should be noted that, for those skilled in the art, several changes and improvements can be made without departing from the scope of the present invention. These all fall within the scope of protection of the present invention.

[0089] The embodiment of the present invention provides a topology optimization method considering stability under the action of design-related loads, referring to Figure 1 As shown, the method specifically includes the following contents:

[0090] Step S1: Initialize the level set function, input a set of radial basis functions and interpolate the level set function to obtain the interpolation coefficients of this set of radial basis functions;

[0091] Step S2: construct a topology optimization model with the structural shape as the design variable, the minimization of the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure as the optimization goal, and the volume of the material as the optimization constraint;

[0092] Step S3: According to the current level set function, search for the current force boundary contour of the structure, and use it to calculate the equivalent nodal force currently acting on the structure;

[0093] Step S4: Calculate the current objective function, shape derivatives of optimization constraints, and Lagrangian operators, solve the velocity field based on the calculation results to evolve and obtain a new level set function, and update the interpolation coefficients before the radial basis function;

[0094] Step S5: Determine the convergence. If the convergence condition is met, output the optimization results under the design-related loads with consideration of flexibility minimization and proceed to the next stage of optimization. Otherwise, return to step S3.

[0095] Step S6: Using the optimization result obtained in step S5 as the initial design, a topology optimization problem is performed with the goal of minimizing the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor and the material volume fraction as the optimization constraint. Based on the current level set function, the current force boundary contour of the structure is searched, and the equivalent nodal force currently acting on the structure is calculated.

[0096] Step S7: Under the current level set function and the load, calculate the inverse of the first four-order eigenvalue buckling factors of the structure;

[0097] Step S8: Calculate the current objective function, shape derivatives of optimization constraints, and Lagrangian operators, solve the velocity field based on the calculation results to evolve and obtain a new level set function, and update the interpolation coefficients before the radial basis function;

[0098] Step S9: Determine the convergence. If the convergence condition is met, output the stability maximization topology optimization result under the design-related load. Otherwise, return to step S6.

[0099] Specifically, the expression of the topology optimization model mentioned in the steps of the invention for minimizing the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure includes:

[0100]

[0101]

[0102]

[0103]

[0104]

[0105]

[0106] Where J(Ω) is the current objective function value; u is the displacement field of the structure, p is the pressure on the Neumann boundary, p0 is the pressure value; v is the test function field; U ad is the displacement field space; J KS (μ i (Ω)) is the current objective function value; ρ is the relative volume density under the imaginary density model; n is the normal vector on the structure contour, φ is the level set function, V(Ω) = ∫ Ω dx represents the volume of the structure, μ represents the inverse of the eigenvalue buckling factor of the structure, μ i (Ω) represents the inverse of the i-th order eigenvalue buckling factor of the structure.

[0107] The characteristic buckling factor λ of the structure is obtained by the following formula:

[0108] ∫ Ω Ae(ψ)e(ψ)dx+λ∫ Ω Ae(u):q(ψ,ψ)dx=0

[0109] Where q(ψ,ψ)=ψ k,i ψ k,j , ψ is the mode of the structure after eigenvalue buckling occurs. Specifically, ψ k,i ψ k,j represents the inner product sum operation between the derivatives of each side of ψ; in order to simplify the optimization problem into a minimization problem, the eigenvalue buckling factor λ is taken as its inverse μ, and the inverse eigenvalue buckling factor μ of the structure is obtained by the following formula:

[0110] ∫ Ω Ae(u):q(ψ,ψ)dx+μ∫ Ω Ae(ψ)e(ψ)dx=0.

[0111] In step S4, in order to smooth the optimization process and the obtained results, the inverse of the first four order buckling factors is subjected to KS aggregation and used as the objective function:

[0112]

[0113] Among them, ρ is a parameter in the KS aggregation equation, which is used to make the result of KS aggregation approach the inverse of the first-order buckling factor μ1. The role of the KS aggregation equation is to ensure the smooth progress of the optimization process.

[0114] In step S4, the evolution of the level set function is driven based on the improved HJ equation. The specific formula is as follows:

[0115]

[0116] in, here Refers to a set of known CSRBF radial basis functions, α(t) represents the The interpolation coefficients, is a regularization term used to control the excessive differentiation of the internal structure during the optimization process, v n is the velocity field, which controls the deformation of the structure during the optimization process.

[0117] Solving the velocity field v based on shape derivative theory n , for the topology optimization problem of minimizing the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure, the Lagrangian equations of each model are listed to solve the adjoint fields of each problem:

[0118]

[0119]

[0120] Among them, the topology optimization problem of minimizing structural flexibility is a self-adjoint problem. For the topology optimization problem of minimizing the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factors of the structure, the adjoint variable χ of the inverse of the eigenvalue buckling factor is obtained by the following formula:

[0121]

[0122] in, Used to normalize the eigenvalue vector. is the displacement field under the current problem; e(χ) is the strain tensor accompanying the field χ.

[0123] The shape derivatives of the Lagrange equations for each model are:

[0124]

[0125]

[0126] Where div(*) represents the divergence of *; θ represents the descending direction of the shape derivative.

[0127] The material volume fraction of the structure is controlled by the augmented Lagrangian operator to reach the optimization constraint level. The augmented Lagrangian operator is obtained by the following formula:

[0128]

[0129] Where V1 is the material volume fraction of the initial structure, μ is a fixed parameter, and V0 is the volume fraction constraint value; n r is the number of iterations of the Lagrangian operator form transformation, Λ k is the Lagrangian operator value at the kth iteration step, γ k Represents the value of γ in the γth iteration step. The role of γ is to control the material volume fraction of the structure to approach the volume fraction constraint value. k The optimization process continues to change as shown by the following formula:

[0130] γ k+1 =min(γ k +Δγ,γ max ),k>n r

[0131] Among them, Δγ is the increment of γ in the optimization process, γ max is the upper limit of γ.

[0132] Velocity field v for each problem n Obtained by the following formula:

[0133]

[0134]

[0135] in,

[0136]

[0137] Where δ(φ) is the Dirac function obtained by the following formula:

[0138]

[0139] Among them, Δ is greater than 5, and Ω1 and Ω2 represent the loaded and unloaded areas respectively;

[0140] An implicit boundary search method is used to determine the load application boundary of the structural design. The search strategy is: perform directional scanning row by row or column by column along the direction of the external load, and perform lateral scanning on each node in the direction perpendicular to the scanning direction in turn until the node φ ≥ 0, and stop scanning. The node is marked as the main mark point, and then the adjacent unit nodes are marked as secondary mark points. Only at the main and secondary mark points, the hydrostatic pressure coefficient p0(x) is not 0; it is used to distinguish the areas with and without load Ω1 and Ω2. The equivalent nodal force of the loaded area is finally determined by the following formula:

[0141]

[0142] Among them, N i (ξ,η) is the shape function, ξ and η are the horizontal and vertical coordinates in the isoparametric coordinate system, and their values range from -1 to 1.

[0143] The present invention also provides a topology optimization system that considers stability under design-related loads. The topology optimization system can be implemented by executing the process steps of the topology optimization method. That is, those skilled in the art can understand the topology optimization method as a preferred embodiment of the topology optimization system. The system specifically includes the following:

[0144] Module M1: Initialize the level set function, input a set of radial basis functions and interpolate the level set function to obtain the interpolation coefficients of this set of radial basis functions;

[0145] Module M2: Construct a topology optimization model with the structural shape as the design variable, the minimization of the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure as the optimization objective, and the volume of the material as the optimization constraint;

[0146] Module M3: Based on the current level set function, search for the current force boundary contour of the structure and use it to calculate the equivalent nodal force currently acting on the structure;

[0147] Module M4: Calculates the current objective function, shape derivatives of optimization constraints, and Lagrangian operators. Based on the calculation results, it solves the velocity field to evolve and obtain a new level set function, and updates the interpolation coefficients before the radial basis function.

[0148] Module M5: Determine the convergence. If the convergence condition is met, the optimization results under the design-related loads with consideration of flexibility minimization are output and the next stage of optimization is carried out. Otherwise, return to module M3.

[0149] Module M6: Using the optimization results obtained in Module M5 as the initial design, the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor is minimized as the goal, and the material volume fraction is used as the optimization constraint for the topology optimization problem. Based on the current level set function, the current force boundary contour of the structure is searched, and the equivalent nodal force currently acting on the structure is calculated.

[0150] Module M7: Calculate the inverse of the first four eigenvalue buckling factors of the structure under the current level set function and load.

[0151] Module M8: Calculates the current objective function, shape derivatives of optimization constraints, and Lagrangian operators. Based on the calculation results, it solves the velocity field to evolve and obtain a new level set function, and updates the interpolation coefficients before the radial basis function.

[0152] Module M9: Determine the convergence. If the convergence conditions are met, output the stability maximization topology optimization results under the design-related loads. Otherwise, return to module M6.

[0153] Next, the present invention will be described in more detail.

[0154] Reference Figure 1 As shown, the topology optimization method provided by the present invention specifically includes the following steps:

[0155] (1) Initialize the level set function, define the structure outline and the inner and outer regions of the structure, and use a set of CSRBF radial basis functions to interpolate the level set function to obtain the interpolation coefficients of this set of radial basis functions.

[0156] Specifically, a set of CSRBF radial basis functions The interpolation of the level set function φ(x,t) under the initial structure can be obtained by the following formula:

[0157]

[0158] Among them, φ(x,t) is the level set function at time t, is the CSRBF radial basis function, and α(t) is the weight coefficient of different radial basis functions in the interpolation process. The formula is:

[0159]

[0160] in, D is the structural design area, d sp is the support radius of the radial basis function, x i Refers to a sample point in the design area, c is a sufficiently small constant, here it is 10 -4 .

[0161] The level set function φ(x,t) is separated into two parts due to the interpolation of the radial basis function. One is a set of radial basis functions that are only related to the spatial variable x. The other is the weight coefficient α(t) that is only related to the time variable t. This method enables the level set function φ(x,t) to complete the decoupling of the time variable t and the space variable x. The sample point positions will not change during the optimization process, that is, the radial basis function group is constant, so only the weight coefficient α(t) related to the time variable t changes with time. Therefore, in the parameterized level set method based on radial basis function interpolation, the only variable is the time variable t, which simplifies the Hamilton-Jacobi equation, the partial differential equation driving the evolution of the level set function, into an ordinary differential equation, as follows:

[0162]

[0163] in,

[0164] Furthermore, in order to avoid excessive differentiation of the internal structure during the stability optimization process, the present invention improves the above formula by adding a regularization term, and the formula is:

[0165]

[0166] in, is the regularization term, and τ is the parameter that controls the degree of internal differentiation of the structure.

[0167] (2) Determine the material parameters and boundary conditions of the structure, and determine the material interpolation model for the stiffness matrix and initial stress matrix of the structure. Specifically, the model formulas are:

[0168] E(ρ)=E min +ρ(E0-E min )

[0169]

[0170] Among them, E0 is the elastic modulus of the solid element material, E min is the elastic modulus of the blank unit material, ρ is the relative volume density under the imaginary density model, and the method for determining the relative volume density under the imaginary density model is as follows: Figure 2 , p m is the penalty factor, and its value range is 0.5 to 1.

[0171] (3) Entering the optimization iteration, the design-related load application boundary of the structure is detected through the implicit expression of the zero level set of the level set function on the structural contour, and the equivalent nodal force in the finite element analysis of the structure is calculated.

[0172] Specifically, the present invention provides a search strategy for the boundaries of structural design-related load application: perform a directional scan row by row or column by column along the direction of the external load, then perform a lateral scan of each node in the scan direction perpendicular to the scan direction until φ ≥ 0. The scan is stopped and marked as a primary marker point. Adjacent unit nodes are then marked as secondary marker points. The hydrostatic pressure coefficient p0(x) is non-zero only at the primary and secondary marker points. The basic concept behind the search for the boundaries of structural design-related load application is illustrated in Figures 3(a) and 3(b).

[0173] In order to accurately and conveniently apply design-related loads to the load-bearing boundaries, the calculation method of equivalent nodal forces in finite element analysis can be obtained:

[0174]

[0175] Where d is the design area, δ ‘ (φ) is the Dirac function, where the Dirac function δ ‘ (φ) is taken as:

[0176]

[0177] Here, ρ ranges from 0.35 to 0.5. The expression of the design-related load at any point in the design area D is:

[0178]

[0179] Furthermore, for the sake of simplicity of calculation, the present invention uses isoparametric elements to perform finite element analysis and calculation. The equivalent nodal force formula at any point in the design area is:

[0180]

[0181] Among them, n i (ξ,η) is the shape function, ξ and η are the horizontal and vertical coordinates in the isoparametric coordinate system, ranging from -1 to 1. The numerical solution is obtained by Gaussian integral:

[0182]

[0183] Where NG is the order of Gaussian integral.

[0184] (4) Based on the parameterized level set method, a topology optimization problem model considering structural stability under the action of design-related loads is established. The displacement field, stress field, and eigenvalue buckling factor proposed in the model are obtained by finite element analysis and eigenvalue buckling analysis calculation methods. Then, based on the shape derivative theory combined with the augmented Lagrangian operator, the velocity field that drives the evolution of the level set function is solved, the level set function is evolved, and the interpolation coefficient before the radial basis function is updated. The specific optimization process is to first perform topology optimization to minimize the structural flexibility. After obtaining the optimal result, the optimal structure at this time is used as the initial structure of the KS aggregation equation value minimization optimization problem of the inverse of the first four-order eigenvalue buckling factor of the structure to obtain the final result.

[0185] Specifically, the model expression of the topology optimization problem considering structural stability under the action of design-related loads can be expressed in stages as follows:

[0186]

[0187]

[0188] V(Ω)=∫ Ω dx≤V0

[0189]

[0190]

[0191] V(Ω)=∫ Ω dx≤V0

[0192] Where u is the displacement field of the structure, p is the pressure on the Neumann boundary, n is the normal vector on the structure contour, Φ is the level set function, V(Ω) = ∫ Ω dx represents the volume of the structure, and μ represents the inverse of the eigenvalue buckling factor of the structure.

[0193] The equilibrium state of the structure based on elastic mechanics can be expressed by the following formula:

[0194]

[0195] Wherein, f is the body force inside the structure. Body force is not considered in the present invention, and here the body force f is taken as 0.

[0196] According to the equilibrium equation and Gaussian divergence theorem, the weak form of the equilibrium equation can be obtained as:

[0197]

[0198] Furthermore, the characteristic buckling factor λ of the structure can be obtained by the following formula:

[0199]

[0200] Where: q(ψ,ψ)=ψ k,i ψ k,j , ψ is the mode of the structure after eigenvalue buckling. To simplify the optimization problem, the eigenvalue buckling factor λ is taken as its inverse μ. The inverse eigenvalue buckling factor μ of the structure can be calculated by the following formula:

[0201]

[0202] Furthermore, in order to optimize the process and smooth the results, the inverse of the first four-order buckling factors is KS aggregated and used as the objective function:

[0203]

[0204] Among them, ρ is a parameter in the KS aggregation equation, which is used to make the result of KS aggregation approach the inverse of the first-order buckling factor μ1.

[0205] Furthermore, the present invention solves the velocity field v based on the shape derivative theory. n , for the topology optimization problem of minimizing the structural flexibility and the KS aggregation equation of the inverse of the first four-order eigenvalue buckling factor of the structure, first list the Lagrangian equations of each model to solve the adjoint fields of each problem:

[0206]

[0207]

[0208] For L(Ω,u,q), The displacement field and its adjoint field are used to find their variation, where L(Ω,u,q) is a self-adjoint problem: u=-q. The following mainly discusses By chain derivation, It can be decomposed into the following formula:

[0209]

[0210] in, It can also be expressed as:

[0211]

[0212] Furthermore, in order to obtain the adjoint field of this problem, due to is the stationary point of the Lagrangian function that satisfies the optimality condition of the minimization problem, so the Lagrangian function is about The variations of are expressed as follows:

[0213]

[0214]

[0215] in, yes The descending direction, for the sake of simplicity, ∫ Ω Ae(ψ i )e(ψ i )dx normalization: ∫ Ω Ae(ψ i )e(ψ i )dx=1,ψ i ′ can be obtained by the following formula:

[0216]

[0217] Furthermore, for the Lagrangian function The variation of can be improved as follows:

[0218]

[0219] After simplification, the adjoint field of this problem is obtained. The formula of the adjoint field in the structure Ω is:

[0220]

[0221] Among them, χ is the adjoint field in the structure Ω, which can be obtained by finite element analysis. Then, for the Lagrangian equation L(Ω,u,q), The shape derivatives can be solved:

[0222]

[0223]

[0224] Among them, H = div(n), the present invention aims at the topology optimization problem under the action of design-related loads, and then the item Can be simplified to:

[0225]

[0226]

[0227] Furthermore, the shape derivative of the above Lagrangian function can be changed to:

[0228]

[0229]

[0230] Furthermore, according to the chain rule, The shape derivative formula is:

[0231]

[0232] Furthermore, the design constraints of the topology optimization problem on the material volume fraction are met, and the Lagrangian function can be expressed as:

[0233]

[0234]

[0235] Where Λ1 is the Lagrangian operator for the material volume fraction constraint. Further, for the above Lagrangian function L1(Ω,u,q), The shape derivative of can be expressed as:

[0236]

[0237]

[0238] The present invention proposes a corresponding augmented Lagrangian operator to control the material volume fraction of the structure to reach an optimized constraint level. The augmented Lagrangian operator is obtained by the following formula:

[0239]

[0240] Where V1 is the material volume fraction of the initial structure, μ is a fixed parameter, and γ k The optimization process continues to change as shown by the following formula:

[0241] γ k+1 =min(γ k +Δγ,γ max ),k>n r

[0242] Among them, Δγ is the increment of γ in the optimization process, γ max is the upper limit of γ.

[0243] Furthermore, the velocity fields of the above problems can be expressed as:

[0244]

[0245]

[0246] in,

[0247]

[0248] The present invention converts the line integral into the volume integral by the following formula to facilitate the implementation of the method:

[0249]

[0250] Among them, in the level set function δ(φ) is the Dirac function and can be obtained by the following formula:

[0251]

[0252] Among them, the value of Δ is greater than 5.

[0253] By converting line integrals into volume integrals, the velocity field can be expressed as:

[0254]

[0255]

[0256] In addition, Ω1 and Ω2 represent the areas with and without load, respectively.

[0257] This paper proposes, for the first time, a topology optimization method based on the level set method that considers structural stability under design-related loads. This method utilizes shape derivative theory and the augmented Lagrangian operator to solve for the velocity field under design-related loads, taking into account structural stability. Furthermore, a regularization term is introduced into the HJ equation to control excessive structural differentiation during the optimization process, avoiding the formation of fine structural branches within the structure. This improves the rationality and manufacturability of the final structure.

[0258] Examples:

[0259] In the design of pressure-resistant structures, the structure is subjected to external hydrostatic pressure. This paper gives a two-dimensional example of the design of a pressure-resistant arch structure subjected to water pressure to illustrate the effectiveness of this method. The design area and boundary conditions of this example are given by Figure 4 As shown in the figure, the design area is a rectangle of (6m×4m). The structure is subject to fixed constraints near the left and right ends, with a constraint width of 0.3m. The hydrostatic pressure is 1. The design area is discretized into a 150×90 finite element mesh (a plane stress model is used in this study) using four-node isoparametric elements with unit width and height. The elastic modulus of the material of the structure is taken as 1, and the Poisson's ratio is set to 0.3. Figure 5(a) shows the optimization results under flexibility minimization, and Figure 5(b) shows the optimal stability result when the final result of flexibility minimization is the initial design.

[0260] The first-order buckling mode of the optimization results for flexibility minimization and the stability optimization result are shown in Figure 6. Compared to the optimization results for flexibility minimization, the optimization results for stability maximization generate multiple holes at the maximum deformation point in the first-order buckling mode of the optimization results for flexibility minimization within the optimized structure, which reduces the deformation degree of the structure at the critical section in the first-order buckling mode, thereby improving the stability performance of the structure. The result graph also shows that the internal boundaries and load-bearing boundaries of the structure obtained by the topology optimization based on the parameterized level set method considering the structural stability under the action of design-related loads are smooth, and holes can be automatically generated within the structure during the optimization process.

[0261] The first-order buckling factor, flexibility and volume fraction change curves of the structure during the stability optimization process are respectively given by Figure 7 , Figure 8 and Figure 9 As shown. Since the present invention adopts the augmented Lagrangian operator, the various indicators shown in the figure gradually stabilize during the oscillation until convergence. From the change diagram of the optimization process of the first-order buckling factor, it can be found that after considering the minimization optimization of the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure as the optimization target, the stability of the structure is significantly improved, and the change of the first-order buckling factor during the optimization process is relatively stable. From the change diagram of the optimization process of flexibility, it can be found that the flexibility gradually increases during the optimization process, that is, the stiffness of the structure becomes smaller and smaller, and the stability of the structure shows a negative correlation with the stiffness. However, compared with the increase in the first-order buckling factor, the decrease in the stiffness of the structure is small.

[0262] The embodiments of the present invention provide a topology optimization method and system that consider stability under the action of design-related loads, which can solve the topology optimization problem that considers stability under the action of design-related loads.

[0263] Those skilled in the art will appreciate that, in addition to implementing the system and its various devices, modules, and units provided by the present invention in purely computer-readable program code, it is entirely possible to implement the same functions of the system and its various devices, modules, and units provided by the present invention in the form of logic gates, switches, application-specific integrated circuits, programmable logic controllers, and embedded microcontrollers by logically programming the method steps. Therefore, the system and its various devices, modules, and units provided by the present invention can be considered a hardware component, and the devices, modules, and units included therein for implementing various functions can also be considered as structures within the hardware component; the devices, modules, and units for implementing various functions can also be considered as both software modules implementing the method and structures within the hardware component.

[0264] The above describes specific embodiments of the present invention. It should be understood that the present invention is not limited to the specific embodiments described above, and those skilled in the art may make various changes or modifications within the scope of the claims, which do not affect the essence of the present invention. The embodiments of this application and the features in the embodiments may be combined with each other in any manner unless there is a conflict.

Claims

1. A topology optimization method considering stability under design-related loads, characterized in that: include: Step S1: Initialize the level set function, input a set of radial basis functions and interpolate the level set function to obtain the interpolation coefficients of this set of radial basis functions; Step S2: construct a topology optimization model with the structural shape as the design variable, the minimization of the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure as the optimization goal, and the volume of the material as the optimization constraint; Step S3: According to the current level set function, search for the current force boundary contour of the structure, and use it to calculate the equivalent nodal force currently acting on the structure; Step S4: Calculate the current objective function, shape derivatives of optimization constraints, and Lagrangian operators, solve the velocity field based on the calculation results to evolve and obtain a new level set function, and update the interpolation coefficients before the radial basis function; Step S5: Determine the convergence. If the convergence condition is met, output the optimization results under the design-related loads with consideration of flexibility minimization and proceed to the next stage of optimization. Otherwise, return to step S3. Step S6: Using the optimization result obtained in step S5 as the initial design, a topology optimization problem is performed with the goal of minimizing the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor and the material volume fraction as the optimization constraint. Based on the current level set function, the current force boundary contour of the structure is searched, and the equivalent nodal force currently acting on the structure is calculated. Step S7: Under the current level set function and the load, calculate the inverse of the first four-order eigenvalue buckling factors of the structure; Step S8: Calculate the current objective function, shape derivatives of optimization constraints, and Lagrangian operators, solve the velocity field based on the calculation results to evolve and obtain a new level set function, and update the interpolation coefficients before the radial basis function; Step S9: judging convergence. If the convergence condition is met, outputting the stability maximization topology optimization result under the design-related loads; otherwise, returning to step S6; In step S2, the expression of the topology optimization model for minimizing the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure includes: Where J(Ω) is the current objective function value; u is the displacement field of the structure, p is the pressure on the Neumann boundary, p0 is the pressure value; v is the test function field; U ad is the displacement field space; J KS (μ i (Ω)) is the current objective function value; ρ is the relative volume density under the imaginary density model; n is the normal vector on the structure contour, Φ is the level set function, V(Ω)=∫ Ω dx represents the volume of the structure, μ represents the inverse of the eigenvalue buckling factor of the structure, μ i (Ω) represents the inverse of the i-th order eigenvalue buckling factor of the structure.

2. The topology optimization method considering stability under design-related loads according to claim 1, characterized in that: The characteristic buckling factor λ of the structure is obtained by the following formula: Where q(ψ,ψ)=ψ k ,iψ k,j , ψ is the mode of the structure after eigenvalue buckling occurs. Specifically, ψ k,i ψ k,j represents the inner product sum operation between the derivatives of each side of ψ; in order to simplify the optimization problem into a minimization problem, the eigenvalue buckling factor λ is taken as its reciprocal μ, and the reciprocal eigenvalue buckling factor μ of the structure is obtained by the following formula: ∫ Ω Ae(u):q(ψ,ψ)dx+μ∫ Ω Ae(ψ)e(ψ)dx=0.

3. The topology optimization method considering stability under design-related loads according to claim 1, characterized in that: In step S4, in order to smooth the optimization process and the obtained results and ensure the smooth progress of the optimization process, the inverse of the first four-order buckling factors is subjected to KS aggregation and used as the objective function: Among them, ρ is a parameter in the KS aggregation equation, which is used to make the result of KS aggregation approach the inverse of the first-order buckling factor μ1.

4. The topology optimization method considering stability under design-related loads according to claim 1, characterized in that: In step S4 and step S8, the evolution of the level set function is driven based on the improved HJ equation, and the specific formula is as follows: in, here Refers to a set of known CSRBF radial basis functions, ɑ(t) represents the The interpolation coefficients, is a regularization term used to control the excessive differentiation of the internal structure during the optimization process, v n is the velocity field, which controls the deformation of the structure during the optimization process.

5. The topology optimization method considering stability under design-related loads according to claim 4, characterized in that: Solving the velocity field v based on shape derivative theory n , for the topology optimization problem of minimizing the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure, the Lagrangian equations of each model are listed to solve the adjoint fields of each problem: Among them, the topology optimization problem of minimizing structural flexibility is a self-adjoint problem. For the topology optimization problem of minimizing the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factors of the structure, the adjoint variable χ of the inverse of the eigenvalue buckling factor is obtained by the following formula: in, Used to normalize the eigenvalue vector, is the displacement field under the current problem; e(χ) is the strain tensor accompanying the field χ.

6. The topology optimization method considering stability under design-related loads according to claim 5, characterized in that: The shape derivatives of the Lagrange equations of each model in the method are: Where div(*) represents the divergence of *; θ represents the descending direction of the shape derivative; The material volume fraction of the structure is controlled by the augmented Lagrangian operator to reach the optimization constraint level. The augmented Lagrangian operator is obtained by the following formula: Where V1 is the material volume fraction of the initial structure, μ is a fixed parameter, and V0 is the volume fraction constraint value; n r is the number of iterations of the Lagrangian operator form transformation, Λ k is the Lagrangian operator value at the kth iteration step, γ k Represents the value of γ at the kth iteration step. The role of γ is to control the material volume fraction of the structure to approach the volume fraction constraint value. k The optimization process continuously changes as shown by the following formula: c k+1 =min(γ k +Δγ,γ max ),k>n r Among them, Δγ is the increment of γ in the optimization process, γ max is the upper limit of γ; The velocity field v of each problem in the method n Obtained by the following formula: in, Where δ(φ) is the Dirac function obtained by the following formula: Among them, Δ is greater than 5, and Ω1 and Ω2 represent the loaded and unloaded areas respectively; An implicit boundary search method is used to determine the load application boundary related to the structural design. A directional scan is performed row by row or column by column along the direction of the external load. Each node is scanned laterally in the direction perpendicular to the scanning direction in turn until φ ≥ 0. The scan is stopped and marked as the main mark point. The adjacent unit nodes are then marked as secondary mark points. Only at the main and secondary mark points is the hydrostatic pressure coefficient p0(x) not 0; this is used to distinguish between the load-free areas Ω1 and Ω2. The equivalent nodal force of the loaded area is finally determined by the following formula: Among them, N i (ξ,η) is the shape function, ξ and η are the horizontal and vertical coordinates in the isoparametric coordinate system, and their values range from -1 to 1.

7. A topology optimization system considering stability under design-related loads, characterized in that: include: Module M1: Initialize the level set function, input a set of radial basis functions and interpolate the level set function to obtain the interpolation coefficients of this set of radial basis functions; Module M2: Construct a topology optimization model with the structural shape as the design variable, the minimization of the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure as the optimization objective, and the volume of the material as the optimization constraint; Module M3: Based on the current level set function, search for the current force boundary contour of the structure and use it to calculate the equivalent nodal force currently acting on the structure; Module M4: Calculates the current objective function, shape derivatives of optimization constraints, and Lagrangian operators. Based on the calculation results, it solves the velocity field to evolve and obtain a new level set function, and updates the interpolation coefficients before the radial basis function. Module M5: Determine the convergence. If the convergence condition is met, the optimization results under the design-related loads with consideration of flexibility minimization are output and the next stage of optimization is carried out. Otherwise, return to module M3. Module M6: Using the optimization results obtained in Module M5 as the initial design, the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor is minimized as the goal, and the material volume fraction is used as the optimization constraint for the topology optimization problem. Based on the current level set function, the current force boundary contour of the structure is searched, and the equivalent nodal force currently acting on the structure is calculated. Module M7: Calculate the inverse of the first four eigenvalue buckling factors of the structure under the current level set function and load. Module M8: Calculates the current objective function, shape derivatives of optimization constraints, and Lagrangian operators. Based on the calculation results, it solves the velocity field to evolve and obtain a new level set function, and updates the interpolation coefficients before the radial basis function. Module M9: Determine the convergence. If the convergence conditions are met, output the stability maximization topology optimization results under the design-related loads. Otherwise, return to module M6. In step S2, the expression of the topology optimization model for minimizing the structural flexibility and the KS aggregation equation value of the inverse of the first four-order eigenvalue buckling factor of the structure includes: Where J(Ω) is the current objective function value; u is the displacement field of the structure, p is the pressure on the Neumann boundary, p0 is the pressure value; v is the test function field; U ad is the displacement field space; J KS (μ i (Ω)) is the current objective function value; ρ is the relative volume density under the imaginary density model; n is the normal vector on the structure contour, Φ is the level set function, V(Ω)=∫ Ω dx represents the volume of the structure, μ represents the inverse of the eigenvalue buckling factor of the structure, μ i (Ω) represents the inverse of the i-th order eigenvalue buckling factor of the structure.

8. A method for preparing a pressure-resistant structure, characterized in that: A pressure-resistant structure is prepared according to a pressure-resistant structure corresponding to a stability-maximizing topology optimization result under design-related loads output by the method of any one of claims 1 to 6.

9. A pressure-resistant structure, characterized in that: The pressure-resistant structure is a pressure-resistant structure prepared by outputting a pressure-resistant structure corresponding to a topology optimization result of maximizing stability under design-related loads according to the method of any one of claims 1 to 6.

Images