An elasto-plastic isogeometric topology optimization method based on BESO

Abstract

The application discloses an elastic-plastic isometric topology optimization method based on BESO, and belongs to the isometric topology optimization field.The method considers isotropic hardening of material through elastic-plastic modeling, equivalently converts stress of the material into Von Mises stress, and allows elastic and plastic deformation of the structure under external loading; the optimization process is closer to actual engineering requirements, and the structure has sufficient strength and stability under complex loading conditions; the isometric analysis technology is adopted, and structure response of the model can be more accurately obtained; the elastic-plastic response is comprehensively considered, and corresponding adjustment is made in the topology optimization, so that the topology layout of the structure is adjusted according to the elastic-plastic properties of the material to ensure the best performance; the BESO optimization method, the isometric analysis and the elastic-plastic material are closely combined, and the vacancy of the previous elastic-plastic isometric topology optimization method is made up, and the plastic behavior of the material can be considered in the optimization.

Description

Technical Field

[0001] This invention belongs to the field of isogeometric topology optimization, and more specifically, relates to an elastoplastic isogeometric topology optimization method based on BESO. Background Technology

[0002] With the development of computer technology, topology optimization has become a powerful tool for achieving innovative designs. Bi-directional Evolutionary Structural Optimization (BESO), a topology optimization method, is an algorithm for finding the optimal material distribution in a structure. Its goal is to find the optimal material distribution in a given space to satisfy a set of performance and constraint conditions, while minimizing material usage through a binary material distribution.

[0003] Isogeometric analysis (IGA) differs from traditional finite element analysis (FEA). It directly uses non-uniform rational B-splines (NURBS) as geometric representations for analysis, eliminating the time-consuming conversion process and avoiding errors introduced during the conversion. This enables seamless integration between computer-aided design (CAD) and computer-aided engineering (CAE).

[0004] The advent of IGA has significantly improved the efficiency and accuracy of numerical simulation and design processes, leading to the widespread development of BESO-based isogeometric topology optimization methods. However, most of these methods revolve around the linear elastic assumption. In practical engineering, numerous material nonlinearities exist, with elastoplastic materials being a common example. When an elastoplastic material deforms under stress to the point where the internal stress exceeds the proportional limit, its stress-strain curve exhibits a nonlinear relationship. As the stress level further increases, irreversible plastic deformation occurs when the stress exceeds the yield stress. Currently, most isogeometric topology optimization methods do not consider the true characteristics of these materials. This limitation makes previous methods difficult to adapt to practical engineering applications, especially when considering the deformation and plastic behavior of structural materials during loading, resulting in optimization results that are significantly inconsistent with reality. Summary of the Invention

[0005] To address the aforementioned deficiencies or improvement needs of existing technologies, this invention provides a BESO-based elastoplastic isogeometric topology optimization method. This solves the problem that existing isogeometric topology optimization methods focus on the rigidity and static loading conditions of the structure, without fully considering the differences between elastic and plastic behaviors, resulting in the final optimization results being seriously inconsistent with the actual situation.

[0006] To achieve the above objectives, according to a first aspect of the present invention, a BESO-based method for elastoplastic and other geometric topology optimization is provided, comprising: S1, the structure to be optimized is divided into multiple grid units by equal geometric meshing, and the density of control points of each unit is used as the design variable to construct the optimization model of the structure to be optimized; S2, determine the target volume of the structure to be optimized obtained in this iteration; S3, apply loads to the structure to be optimized, perform elastoplastic and other geometric analyses on each element to obtain the stiffness matrix and equivalent internal forces of each element, and according to... Calculate the displacement increment of each unit control point ;judge If the value is less than the threshold, proceed to step S5; otherwise, proceed to step S6. The displacements are cumulatively added to the control points of each element, and then elastic-plastic geometric analyses are performed on each element again to update the stiffness matrix and equivalent internal forces. The displacement increments of the control points of each element are then calculated again. until If the displacement is less than the threshold, the displacement increment at this time is added to the displacement of each unit control point; S4. Determine whether the displacement of each unit control point has reached the specified displacement. If yes, select a density threshold based on the sensitivity field of all current unit control points and the corresponding volume constraint value to truncate the current sensitivity field to update the design variables, and then update the structural topology. If no, return to S3. S5, determine whether the volume of the current structural topology has reached the target volume and whether the objective function corresponding to the optimization model has converged. If yes, the topology optimization ends and the optimal structural topology is obtained; otherwise, return to S2.

[0007] According to a second aspect of the present invention, a BESO-based elastoplastic geometric topology optimization system is provided, comprising: a computer-readable storage medium and a processor; The computer-readable storage medium is used to store executable instructions; The processor is configured to read executable instructions stored in the computer-readable storage medium and execute the method as described in the first aspect.

[0008] According to a third aspect of the invention, a computer-readable storage medium is provided, the computer-readable storage medium storing computer instructions for causing a processor to perform the method as described in the first aspect.

[0009] In summary, compared with the prior art, the above-described technical solutions conceived by this invention can achieve the following beneficial effects: Based on the BESO topology optimization framework, this paper utilizes an isotropic hardening model and related flow criteria for elastoplastic modeling. Through stress prediction correction schemes and incremental equations, elastoplastic isogeometric analyses are performed on the design model to obtain the structural response, thereby enabling the optimization design of elastoplastic structures. The design of a BESO-based elastoplastic isogeometric topology optimization method is of great significance for nonlinear isogeometric topology optimization of materials.

[0010] 1. The method provided by this invention, through elastoplastic modeling, considers the isotropic hardening of materials and equates material stress to Von Mises stress, allowing the structure to undergo elastic and plastic deformation under external loading. This makes the optimization process closer to actual engineering needs, ensuring that the structure has sufficient strength and stability under complex loading conditions. The use of isogeometric analysis techniques enables more accurate determination of the model's structural response. By comprehensively considering the elastoplastic response and making corresponding adjustments in topology optimization, the topological layout of the structure will be adjusted according to the elastoplastic properties of the material to ensure optimal performance. The close integration of the BESO optimization method, isogeometric analysis, and elastoplastic materials fills the gap in previous methods for elastoplastic isogeometric topology optimization, allowing the plastic behavior of the material to be considered during optimization. This method has broad application potential for various engineering applications, especially in fields requiring high plasticity and strength.

[0011] 2. The method provided by this invention employs a sensitivity numerical processing approach: two filters are used to smooth the sensitivity, avoiding grid dependence and checkerboard topology phenomena. By averaging, the sensitivity information of each step includes all the information from previous iterations, thus improving the stability of the results. Attached Figure Description

[0012] Figure 1 A flowchart of the BESO-based elastoplastic geometric topology optimization method provided in this embodiment of the invention; Figure 2 This is a schematic diagram of the control points of each element in the isogeometric model mesh topology optimization structure provided in the embodiments of the present invention; Figure 3 In the figures (a) and (b), respectively, it is a schematic diagram of the estimation of the objective function value provided by the embodiment of the present invention and a schematic diagram of the influence of displacement control on the objective function; Figure 4 This is a flowchart of solving the constitutive matrix in elastoplastic and other geometric analyses provided in this embodiment of the invention; Figure 5 A schematic diagram of isotropic hardening of two-dimensional units provided in an embodiment of the present invention; Figure 6 This is a schematic diagram illustrating density filtration and NURBS filtration provided in an embodiment of the present invention. Detailed Implementation

[0013] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention. Furthermore, the technical features involved in the various embodiments of this invention described below can be combined with each other as long as they do not conflict with each other.

[0014] This invention provides a BESO-based method for optimizing elastoplastic and other geometric topologies, such as... Figure 1 As shown, it includes: S1. Divide the structure to be optimized into multiple grid cells by equal geometric meshing. Use the density of control points of each cell as the design variable to construct an optimization model for the structure to be optimized.

[0015] The model to be designed is geometrically meshed according to the initially set mesh size to obtain a geometrically meshed model. The geometrically meshed model is then divided to obtain the coordinates and numbers of the corresponding control points and Gaussian points in each element.

[0016] Preferably, the refinement method used in this invention is H-refinement, which refines the design domain according to the initially set mesh size to obtain a uniform geometric mesh model of the specified mesh size, and obtains the coordinates and numbers of the corresponding control points and Gaussian points in each element based on the characteristics of the uniform geometric mesh model, such as... Figure 2 As shown.

[0017] The design objective of this method is to minimize the structural compliance, which is defined by the integral of the external force over a given specified displacement range. This is equivalent to maximizing the mechanical work consumed by the model during deformation. The material volume fraction of the design domain is also specified.

[0018] The mathematical expression for the optimization model of the structure to be optimized is:

[0019] in, , , , and Let the external force vector be an unknown field displacement, the displacement vector be a specified vector, the unbalanced force vector be an unknown vector, and the material volume fraction be a specified vector. For the target value size, For unit density, nel This represents the number of units.

[0020] The convergence condition, filtering radius, and evolution ratio of the objective function are determined based on actual needs.

[0021] The convergence condition can be expressed by the following relation:

[0022] in, and N These represent the allowable convergence error and a positive integer, respectively; in this invention, they are taken as 0.001 and 5, respectively. i This represents the current iteration number.

[0023] S2, determine the target volume of the structure to be optimized obtained in this iteration.

[0024] The target volume for this optimization is determined based on the evolution ratio and the current design domain volume. The target volume for this optimization is determined based on the evolution ratio and the current design domain volume. The current optimization iteration (the...) j The volume fraction of (times) is calculated as follows:

[0025] In the formula, er Indicates the evolutionary ratio. This represents the desired target volume. It's important to note that the actual volume obtained during the iteration process usually deviates slightly from the specified target value, but once the desired material volume fraction is reached, the optimization algorithm will only change the topology, and the volume fraction will remain constant.

[0026] S3, apply loads to the structure to be optimized, and perform elastoplastic and other geometric analyses on each element to obtain the stiffness matrix of each element and the equivalent internal forces at the control points of each element.

[0027] Based on the model, the structural response should be modeled and analyzed elastoplastically for each element to obtain the elastoplastic constitutive matrix for each element and the equivalent Von Mises stress at each Gaussian point. The equivalent internal forces at the control points and the element stiffness matrix of each element are obtained by using the Gaussian point number and the mapping relationship between elements and control points.

[0028] In step S3, the load is applied in a segmented manner. In each loading step, the elastoplastic analysis and the solution of internal forces and element stiffness matrices are performed in the following manner: 1) Obtain the displacements at all control points in each element based on the structural response, and solve for the element strain by combining the partial derivatives of the shape functions at all control points in each element.

[0029] 2) Substitute the element strain into the plastic solver and solve the elastic-plastic matrix of the element and the Von Mises stress at the Gauss point using the prediction correction method.

[0030] 3) Based on the solved elastoplastic matrix and Von Mises stress, solve for the equivalent internal forces at the control points of the corresponding element and the element stiffness matrix.

[0031] Specifically, each element is a two-dimensional element, and the strain of each two-dimensional element is solved according to the following relationship:

[0032] in, It is a shape function in two-dimensional physical coordinate space. n Represents the total number of control points in a building unit. and Under two-dimensional conditions u and v Displacement at different control points in two directions Here is the displacement matrix of the element control points. Here is the strain matrix.

[0033] Isogeometric analysis was performed using segmented load application and Newton's iteration method, such as... Figure 3 As shown in (a) and (b) in the figure; the hardening model of the material is an isotropic hardening model, as follows. Figure 5 As shown, the relevant yield criterion and flow criterion can be expressed by the following formula:

[0034] Where, σ Y , κ and e p Let A be the yield strength, the hardening parameter for yield surface extension, and the effective plastic strain, respectively. Let A be the strain hardening function (the local slope of the uniaxial stress / plastic strain curve).

[0035] Perform elastoplastic analysis and solve for the elastoplastic constitutive matrix. And element stress, the specific process is as follows Figure 4 .

[0036] For a two-dimensional plane problem, strain and its components can be expressed as:

[0037] in, , , .

[0038] Assuming the strain is entirely elastic, the theoretical magnitude of the experimental stress is:

[0039] in, Let be the constitutive matrix under elastic conditions. Let be the initial stress magnitude used in this calculation. Ignoring kinematic hardening, the theoretical experimental shear stress can be expressed as:

[0040] Using the von Mises yield criterion, the final equivalent stress value is obtained, calculated as follows:

[0041] The yield function can be expressed in the following form:

[0042] in, The initial yield strength, a It is a flow vector.

[0043] like f <0 indicates that the element state is elastic, and the element stress, plastic strain, and stiffness matrix are expressed as follows:

[0044] like f A value greater than 0 indicates that the element is in a plastic state, in which case the experimental stress and experimental shear stress are recalculated. The state prior to strain in this iteration is selected as the reference value, and the experimental stress, experimental shear stress, and element stress are calculated as follows:

[0045]

[0046]

[0047] in, For the recalculated experimental stress, For the recalculated shear stress, This represents the updated element stress magnitude.

[0048] After the calculation is completed, the final calculated element stress value is compared with the initial yield strength, and calculations are performed based on their relative magnitudes. r The magnitude of the value. If , ,like ,

[0049] The intermediate variable symbols involved in this step have no actual physical meaning; they serve only as a medium for intermediate derivation. Some variables are defined as follows:

[0050]

[0051] The von Mises equivalent stress can be obtained through the following relationship:

[0052]

[0053]

[0054]

[0055]

[0056]

[0057]

[0058] Stress magnitude calculation:

[0059] Plastic strain calculation:

[0060] To ensure that the final yield point is on the yield surface, a scaling factor is added. Adjust the stress level:

[0061]

[0062] Where D is the elastic constitutive matrix. It is an elasto-plastic matrix, which changes during the iteration process. H The strain hardening function is... Vector Let A be the flow vector. A is the local slope of the uniaxial stress / plastic strain curve, which can be expressed by the equation... Confirmed, among which It is the hardening modulus, also known as the tangential modulus. E It is the elastic modulus.

[0063] The equivalent internal forces and element stiffness matrices at the control points are solved according to the following relationships:

[0064]

[0065] Where N is the number of Gaussian points in the cell. This represents the weighting factor corresponding to the Gaussian point. Let Von Mises stress be the stress at the Gauss point. It is an elastic-plastic matrix. For the physical domain, For the integration domain, These represent the transformation relationships from the NURBS parameter space to the physical space and from the integral parameter space to the NURBS parameter space, respectively. To further optimize the elastoplastic matrix, the transformation relationship is... This can be expressed by the following relation:

[0066]

[0067] in, For physical space coordinates, For parameter space coordinates, These are the coordinates in the integration space.

[0068] S3, apply loads to the structure to be optimized, perform elastoplastic and other geometric analyses on each element to obtain the stiffness matrix and equivalent internal forces of each element, and according to... Calculate the displacement increment of each unit control point ;judge If the value is less than the threshold, proceed to step S5; otherwise, proceed to step S6. The displacements are cumulatively added to the control points of each element, and then elastic-plastic geometric analyses are performed on each element again to update the stiffness matrix and equivalent internal forces. The displacement increments of the control points of each element are then calculated again. until If the displacement is less than the threshold, the displacement increment at this point is added to the displacement of each element control point; where K is the overall stiffness matrix obtained by coupling the stiffness matrices of each element. The overall internal and external force residuals are obtained by coupling the internal and external force residuals of each element; Specifically, the element stiffness matrices of each element are assembled to obtain the overall stiffness matrix, and the structural response of this iteration is solved using equilibrium equations until the internal and external force residual vectors are obtained. .

[0069] That is, assemble the overall stiffness matrix K according to the degree of freedom number of the control points, and then assemble it according to the equilibrium equations. Solve for the displacement increment and determine whether the residual R of the internal and external force vectors satisfies the condition. (Under normal circumstances) If the condition is met, proceed to step S5; otherwise, accumulate the displacement increment obtained in this iteration into the control point displacement, perform elastoplastic analysis again to update the equivalent internal forces at each element control point, solve for the displacement increment again, and determine whether the condition is met again. until satisfied .

[0070] It is worth noting that when performing the plastic and other geometric analyses in S3 for the first time, the equivalent internal forces of each element are zero.

[0071] Solve the equilibrium equations to obtain the structural response. The specific steps are as follows: Calculate the element stiffness matrix And assembled into the overall stiffness matrix, the element stiffness matrix is ​​calculated as follows:

[0072] The element stiffness matrices are assembled into the global stiffness matrix K based on the degree of freedom number.

[0073] The degrees of freedom are divided into boundary degrees of freedom (index D, related to Dirichlet boundary conditions) and residual degrees of freedom (index F), with displacement increments. The stiffness matrix can be expressed as:

[0074] Calculate the displacement of the remaining degrees of freedom This is used as the displacement increment.

[0075] If this is the first iteration, the displacement increment can be calculated using the following formula:

[0076] If this is not the first iteration, then: ,

[0077] in, j This represents the number of iterations. Finally, the corrected displacement increment is added to the displacement vector:

[0078] Calculate each integration point m The strain increment and stress increment are calculated, and the stress increment at each integration point is added to... :

[0079]

[0080] The equivalent internal forces are then calculated, as shown below:

[0081] Finally, compare whether the residuals R of the internal and external force vectors satisfy the condition. Under normal circumstances .

[0082] S4. Determine whether the displacement of each unit control point has reached the specified displacement. If yes, select a density threshold based on the sensitivity field of all current unit control points and the corresponding volume constraint value to truncate the current sensitivity field to update the design variables, and then update the structural topology configuration, and proceed to S5. If no, return to S3 to apply the load again so that the displacement of each unit control point gradually reaches the specified displacement.

[0083] Solve for the sensitivity value and determine whether the load on the specified degree of freedom reaches the specified load (determined by judging whether the displacement of the controlled degree of freedom reaches the specified displacement). If it does, filter and average the sensitivity value to establish a global control point sensitivity field. Use the threshold parameters for material addition and removal obtained by the OC iterative algorithm. and Update the topology variables and proceed to the next step; otherwise, increase the load and return to S3 to continue the elastoplastic analysis.

[0084] The formula for calculating the sensitivity of each unit is as follows:

[0085] in, and For the introduced Lagrange multipliers, and , , , The element strain matrix, For von Mises stress, For the physical domain corresponding to the model, This is the IGA cell domain.

[0086] The obtained sensitivity values ​​are filtered and averaged to establish a global control point sensitivity field. The relevant filtering and truncation operations are as follows: Figure 6 As shown, topology optimization is then performed based on the topology sensitivity values. This step mainly includes the following two parts: (1) Filter and average the sensitivity values ​​to obtain the global control point sensitivity field.

[0087] First, sensitivity filtering is applied to the unit sensitivity values ​​to avoid checkerboard patterns and grid dependency issues, calculated according to the following formula:

[0088] in, nelNumber of units, linear weights Calculate according to the following formula:

[0089] in, and These are the radius of the sensitivity filter and the unit. e and k The distance between centers.

[0090] Then, a NURBS filter is used to filter the sensitivity, obtaining the sensitivity value based on the control point. The specific relationship is as follows:

[0091] in, , ei and eij They represent the first i The variables of the first control point, subject to the first i The set of units affected by each control point and ei The j Units. These are NURBS basis functions.

[0092] Finally, the control point sensitivity values ​​are averaged to ensure that the updated sensitivity values ​​include all sensitivity information from previous iterations. The averaging process is performed according to the following formula:

[0093] in, j This represents the number of optimizations performed.

[0094] At this point, the global control point sensitivity field has been established.

[0095] (2) Obtain the threshold parameters for material addition and removal using the OC iterative algorithm. and To update the topological variables, the iterative algorithm takes the following steps: ① Order , It can be easily obtained from Confirm. For example, if there are 1000 elements in the design domain, and That is, when the target volume is 0.725, For a design with 725 elements of density 1, the density of the other elements should be 0. .

[0096] ② Calculate the allowable volume ratio ( AR This is defined as the number of added cells divided by the total number of cells in the current design. If (in er If the specified maximum volume addition ratio is true, then skip step ③. Otherwise, recalculate according to step ③. and .

[0097] ③ The sensitivity values ​​of empty cells (cells with a density of 0) are first sorted to calculate the sensitivity values. The number of units switching from 0 to 1 will be equal to er Multiply by the total number of units in the current design. It's the sensitivity of the element immediately following the last added unit. Then determine... , such that the volume removed is equal to ( +Add the volume of the unit). For example, if there are 1000 elements in the design, and , er If the value is 0.02, then the number of units switching from 0 to 1 is 1000. 0.02 = 20. Sort the sensitivity values ​​of empty cells from largest to smallest. Based on this order, obtain the sensitivity of the element immediately following the 20th added cell. This determines the sensitivity of the cell. After adding the cells, calculate ( (+volume of added cells) 1000 gives the number of units that need to be removed. l Sort the real unit sensitivities from smallest to largest to obtain the first... l The sensitivity of each element is used to determine... .

[0098] The topology variables are updated according to the following scheme:

[0099] This scheme shows that when the sensitivity of the physical control point is less than Then delete the control point and restore the sensitivity to be greater than 1. The void control points. After the OC iterative update algorithm converges, the density of all element control points is obtained, and the density of each element is calculated by interpolation.

[0100] S5, determine whether the volume of the current structural topology has reached the target volume and whether the objective function corresponding to the optimization model has converged. If yes, the topology optimization ends and the optimal structural topology is obtained; otherwise, return to S2 to enter the next iteration.

[0101] Determine whether the material has reached the target volume and whether the target values ​​of the adjacent iterations meet the convergence condition (i.e., whether the objective function of the optimization model has converged iteratively). If so, stop the iteration and enter the visualization program; otherwise, return to S2 to enter the next iteration.

[0102] This invention provides a BESO-based elastoplastic geometric topology optimization system, comprising: a computer-readable storage medium and a processor; The computer-readable storage medium is used to store executable instructions; The processor is configured to read executable instructions stored in the computer-readable storage medium and execute the method as described in any of the above embodiments.

[0103] This invention provides a computer-readable storage medium storing computer instructions that cause a processor to perform the method described in any of the above embodiments.

[0104] Those skilled in the art will readily understand that the above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A BESO-based method for elastoplastic isogeometric topology optimization, characterized in that, include: S1, the structure to be optimized is divided into multiple grid units by equal geometric meshing, and the density of control points of each unit is used as the design variable to construct the optimization model of the structure to be optimized; S2, determine the target volume of the structure to be optimized obtained in this iteration; S3, apply loads to the structure to be optimized, perform elastoplastic and other geometric analyses on each element to obtain the stiffness matrix and equivalent internal forces of each element, and calculate the displacement increment Δu of the control points of each element according to KΔu=R; determine ||R|| ∞ If the displacement is less than the threshold, proceed to S5; otherwise, accumulate Δu to the displacement of each element control point, perform elastoplastic and other geometric analyses on each element again to update the stiffness matrix and equivalent internal forces, and recalculate the displacement increment Δu of each element control point until ||R|| ∞ If the displacement is less than the threshold, the displacement increment at this time is added to the displacement of each unit control point; where K is the overall stiffness matrix and R is the overall internal and external force residual. S4. Determine whether the displacement of each unit control point has reached the specified displacement. If yes, select a density threshold based on the sensitivity field of all current unit control points and the corresponding volume constraint value to truncate the current sensitivity field to update the design variables, and then update the structural topology. If no, return to S3. S5, determine whether the volume of the current structural topology has reached the target volume and whether the objective function corresponding to the optimization model has converged. If yes, the topology optimization ends and the optimal structural topology is obtained; otherwise, return to S2.

2. The method as described in claim 1, characterized in that, The mathematical expression for the optimization model of the structure to be optimized is: Find:ρ=(ρ1,ρ2,...,ρ nel ) T Minimize: subject to:R=F ext -F int ≈0, Where, ρ e ∈{0;1},e=1,2,...,nel,ρ e is the cell density, nel is the number of cells, F ext F is the external force vector. int For equivalent internal forces, u, V req These represent the unknown field displacement, the specified displacement vector, and the specified material volume fraction, respectively, W. C V is the target value. req Let V be the desired target volume.

3. The method as described in claim 1, characterized in that, The target volume of the structure to be optimized is: In j =max(V req ,V j-1 (1-er)) Among them, V req Let j be the desired target volume, j be the number of iterations, and er be the evolution ratio.

4. The method as described in claim 1, characterized in that, In step S3, the step of performing elastoplastic and other geometric analyses on each element to obtain the stiffness matrix of each element and the equivalent internal forces at the control points of each element includes: Solve for the strain of each element by combining the displacements at all control points of each element with the partial derivatives of the shape functions at all control points of each element. Substitute the strain of each element into the plastic solver, and solve the elastic-plastic matrix of each element and the Von Mises stress at the Gaussian point by the prediction correction method; Based on the solved Von Mises stress, the equivalent internal forces at the control points of each element and the element stiffness matrix of each element are determined.

5. The method as described in claim 1 or 4, characterized in that, The equivalent internal forces at each unit control point are: The element stiffness matrix of each element is as follows: Among them, F int For equivalent internal forces, N is the number of Gaussian points in the element, ω is the weighting factor corresponding to the Gaussian point, B is the element strain matrix, σ is the Von Mises stress at the Gaussian point, and J1 and J2 represent the transformation relationships from NURBS parameter space to physical space and from integral parameter space to NURBS parameter space, respectively; D ep Let Ω be the elastoplastic constitutive matrix, and Ω be the physical domain. This is the integration domain.

6. The method as described in claim 1, characterized in that, The sensitivity fields of all unit control points are obtained as follows: After sensitivity filtering of each unit, the sensitivity of each unit center point is used as the sensitivity of each unit center point, and then used as a sample point to fit a NURBS sensitivity field to obtain the initial sensitivity field of all unit control points. According to the relation The initial sensitivity fields of all unit control points are averaged to obtain the sensitivity fields of all unit control points; where j is the iteration number, i is the control point number, and j≥2.

7. A BESO-based elasto-plastic isogeometric topology optimization system, characterized in that, include: Computer-readable storage media and processors; The computer-readable storage medium is used to store executable instructions; The processor is configured to read executable instructions stored in the computer-readable storage medium and execute the method as described in any one of claims 1-6.

8. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores computer instructions for causing a processor to perform the method as described in any one of claims 1-6.

Images