A Multi-Material Topology Optimization Method Based on Univariate Discrete Interpolation Scheme

Abstract

This invention relates to a multi-material topology optimization method based on a single-variable discrete interpolation scheme, belonging to the field of multi-material topology optimization technology. It solves the problems of numerous design variables, low optimization efficiency, unreasonable design results, and difficult-to-manufacture structures in existing multi-material topology optimization techniques. The method includes the following steps: discretizing the design domain using finite element methods and applying boundary conditions and external loads; constructing a multi-material topology optimization model based on a single-variable discrete combined interpolation method, and then obtaining the sensitivity of the multi-material topology optimization model; updating the design variables based on the multi-material topology optimization model and the sensitivity to obtain the optimal topological density field. This achieves efficient and reasonable multi-material topology optimization.

Description

Technical Field

[0001] This invention relates to the field of structural optimization technology, and in particular to a multi-material topology optimization method based on a single-variable discrete interpolation scheme. Background Technology

[0002] Topology optimization refers to the process of calculating the optimal distribution of materials within a given design domain based on certain constraints (such as material volume and mass) to obtain an efficient structural shape, which can effectively reduce structural weight and improve structural performance.

[0003] In recent years, with the rapid development of additive manufacturing technologies (such as 3D printing), structural topology optimization design has become a hot topic in engineering and scientific research. Over the past two or three decades, topology optimization methods have become increasingly diverse, encompassing methods such as Solid Isotropic Material Penalty (SIMP), level set methods, bidirectional evolutionary structural optimization methods, moving deformable void methods, and phase-field methods. Notably, among density-based methods, the SIMP method is highly favored due to its efficiency and ease of implementation. With the continuous advancement of materials science and manufacturing technology, traditional single-material topology optimization methods can no longer meet the demands for multifunctional and composite performance of structures, and topology optimization has gradually expanded into the multi-material domain. Compared to single-material structures, multi-material structures can better meet complex functional requirements and significantly improve structural performance through the rational combination and distribution of materials. However, topology optimization of multi-material structures still faces many challenges, particularly in terms of the number of design variables, the clarity of material interfaces, and the rationality of material allocation.

[0004] In summary, traditional multi-material topology optimization methods suffer from a "combinatorial explosion" problem as the number of design variables increases with the number of candidate materials. The number of design variables in each design unit is proportional to the number of material phases, resulting in numerous design variables and low optimization efficiency. Existing multi-material topology optimization step interpolation models may exhibit unreasonable designs such as material envelopes, which significantly reduces the optimization design space. The feasible region of the design is significantly smaller than that of discrete material schemes, leading to unreasonable design results and making the designed structures difficult to process and manufacture. Summary of the Invention

[0005] Based on the above analysis, the embodiments of the present invention aim to provide a multi-material topology optimization method based on single-variable discrete interpolation, in order to solve the problems of numerous design variables, low optimization efficiency, unreasonable design results, and difficulty in manufacturing the designed structures in existing multi-material topology optimization methods.

[0006] This invention provides a multi-material topology optimization method based on a single-variable discrete interpolation scheme, comprising the following steps:

[0007] The design domain is discretized using finite element methods, and boundary conditions and external loads are applied.

[0008] Based on the univariate discrete combination interpolation method, a multi-material topology optimization model is constructed, and then the sensitivity of the multi-material topology optimization model is obtained;

[0009] Based on the multi-material topology optimization model and sensitivity, the design variables are updated to obtain the optimal topology density field.

[0010] Furthermore, based on the univariate discrete combination interpolation method, a multi-material topology optimization model is constructed, including:

[0011] Based on univariate characteristic functions, design variables are mapped to the topological density field of finite element elements;

[0012] Based on the topological density field of the finite element, the filtered topological density field of the finite element is obtained, and then the normalized topological density field of the finite element is obtained.

[0013] Based on the normalized topological density field of the finite element, the elastic modulus of the candidate material is discretized and interpolated to obtain the elastic modulus of the finite element, and then the elastic tensor of the finite element is obtained.

[0014] A multi-material topology optimization model is constructed based on the elastic tensor of finite element elements.

[0015] Furthermore, the univariate characteristic function is expressed as:

[0016]

[0017] In the formula, This indicates that when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... Topological density of phase candidate materials, This indicates that when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... The Heaviside mapping function for the candidate material, where η represents the penalty function. Θ Indicates the phase number of the candidate material.

[0018] Furthermore, when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... Heaviside mapping function of phase candidate materials Represented as:

[0019]

[0020] In the formula, Indicates the first The threshold parameter for the candidate material, β represents the projection parameter, and tanh() represents the hyperbolic tangent function.

[0021] Furthermore, the elastic modulus and elastic tensor of the finite element are expressed as follows:

[0022]

[0023] In the formula, E[χ(ξ)] and D[χ(ξ)] represent the elastic modulus and elastic tensor of the ξ-th finite element when the design variable is χ(ξ), respectively; Let these represent the elastic modulus and elastic tensor of the empty material, respectively. They represent the first The elastic modulus and elastic tensor of the candidate materials, It indicates the ξ-th finite element. The weighting coefficients of candidate materials.

[0024] Furthermore, the ξ-th finite element element... Weighting coefficients of candidate materials Represented as:

[0025]

[0026] In the formula, These represent the design variables of the ξ-th finite element element as χ(ξ). Normalized topological density of k-phase candidate materials.

[0027] Furthermore, based on the multi-material topology optimization model and sensitivity, the gradient optimization algorithm MMA is used to update the design variables until the multi-material topology optimization model converges, thereby obtaining the optimal design variables and the optimal topology density field.

[0028] Furthermore, when updating design variables using the gradient optimization algorithm MMA, the optimization parameters in the iteration are updated based on an adaptive parameter adjustment strategy; the optimization parameters include penalty parameters and projection parameters.

[0029] Furthermore, the adaptive parameter adjustment strategy includes:

[0030] Before the gradient optimization algorithm MMA iteration, the optimization parameters are initialized, and the counter value is initialized to 0;

[0031] During the MMA gradient optimization algorithm iteration, it is determined whether the projected parameters in the current iteration are less than the set projection parameter threshold. If so, the design variable in the current iteration is taken as the optimal design variable; otherwise:

[0032] Determine if the counter value is greater than 0 in the current iteration.

[0033] If yes, then determine whether the change in the objective function in the current iteration is less than the predetermined value of the change. If yes, then set the counter value to 1 and do not update the optimization parameters in the current iteration; otherwise, do not update the optimization parameters in the current iteration.

[0034] If not, the counter value is incremented by 1, and it is determined whether the current counter value is the set number of times threshold. If yes, the optimization parameters in the current iteration are updated, and the counter value is set to 1; otherwise, the optimization parameters in the current iteration are not updated.

[0035] Furthermore, the projection parameters in the optimized parameters are updated using the following formula:

[0036] β=θ×β′

[0037] In the formula, β′ is the projection parameter in the current iteration, β is the updated projection parameter, and θ represents the first adjustment coefficient.

[0038] The penalty parameter in the optimization parameters is updated using the following formula:

[0039]

[0040] In the formula, η′ is the penalty parameter in the current iteration, η is the penalty parameter after the update, and σ represents the second adjustment coefficient.

[0041] Compared with the prior art, the present invention can achieve at least one of the following beneficial effects:

[0042] This invention provides a multi-material topology optimization method based on a single-variable discrete interpolation scheme. By discretizing the design domain using finite element methods and applying boundary conditions and external loads, a multi-material topology optimization model is constructed using a single-variable discrete combination interpolation method. The sensitivity of the multi-material topology optimization model is then obtained. Based on the multi-material topology optimization model and the sensitivity, the design variables are updated to obtain the optimal topological density field. This solves the problems of numerous design variables, low optimization efficiency, unreasonable design results, and difficulty in manufacturing the designed structure. By constructing a single-variable characteristic function and using a single variable to select multiple discrete physical quantities, the composite function of the reduced-order multi-characteristic function of the single variable is combined with the corresponding material properties for material interpolation. This eliminates the material envelope phenomenon present in existing step-type interpolation models and overcomes the disadvantage of a limited design domain. The number of design variables is independent of the number of candidate material phases, significantly reducing the number of design variables. While achieving similar or even lower objective function values, fewer iterations are needed, improving optimization efficiency. It can achieve a reasonable distribution of materials in the structure while ensuring material mass or volume fraction constraints.

[0043] In this invention, the above-described technical solutions can be combined with each other to achieve more preferred combinations. Other features and advantages of this invention will be set forth in the following description, and some advantages may become apparent from the description or be learned by practicing the invention. The objects and other advantages of this invention can be realized and obtained from what is particularly pointed out in the description and drawings. Attached Figure Description

[0044] The accompanying drawings are for illustrative purposes only and are not intended to limit the invention. Throughout the drawings, the same reference numerals denote the same parts.

[0045] Figure 1 This is a flowchart illustrating the multi-material topology optimization method based on a single-variable discrete interpolation scheme provided in Embodiment 1 of the present invention.

[0046] Figure 2 This is a schematic diagram of the specific process of the multi-material topology optimization method based on a single-variable discrete interpolation scheme provided in Embodiment 1 of the present invention;

[0047] Figure 3 This is a schematic diagram of the conditions for a two-dimensional MBB beam under scenario 1 provided in Embodiment 2 of the present invention;

[0048] Figure 4 This is a schematic diagram of the iterative convergence process of the objective function value and volume fraction of the MBB beam structure under scenario 1 provided in Embodiment 2 of the present invention.

[0049] Figure 5(a) is a comparison of the two-dimensional MBB beam topology optimization results obtained by the proposed method with the DMO method under scenario 1 provided in Embodiment 2 of the present invention when θ = 1.24 and σ = 1.065.

[0050] Figure 5(b) is a comparison of the two-dimensional MBB beam topology optimization results obtained by the proposed method with the DMO method under scenario 1 provided in Embodiment 2 of the present invention when θ = 1.1445 and σ = 1.0247.

[0051] Figure 6 This is a schematic diagram of the conditions for a cantilever beam structure under scenario 2 provided in Embodiment 2 of the present invention;

[0052] Figure 7 This is a graph showing the trend of the objective function and volume fraction as a function of the number of iterations under scenario 2 of embodiment 2 of the present invention, with respect to the quantities of the five materials.

[0053] Figure 8(a) is a comparison between the two-dimensional MBB beam topology optimization results obtained by the proposed method under scenario 2 with five materials provided in Embodiment 2 of the present invention and the DMO method;

[0054] Figure 8(b) is a comparison of the two-dimensional MBB beam topology optimization results obtained by the proposed method under scenario 1 with three materials provided in Embodiment 2 of the present invention and the DMO method;

[0055] Figure 9 This is a schematic diagram of the bridge conditions under scenario 3 provided in Embodiment 2 of the present invention;

[0056] Figure 10 This is a schematic diagram showing the objective function value and the iterative convergence of the volume fractions of the three materials in the multi-material topology optimization process of the bridge structure under scenario 3 provided in Embodiment 2 of the present invention.

[0057] Figure 11 The material distribution results of the three-dimensional bridge structure under scenario 3 provided in Embodiment 2 of the present invention after multi-material topology optimization;

[0058] Figure 12 This is a comparison diagram of the 3D bridge structure topology optimization results and the DMO method under scenario 3 provided in Embodiment 2 of the present invention. Detailed Implementation

[0059] Preferred embodiments of the present invention will now be described in detail with reference to the accompanying drawings, which form part of this application and are used together with the embodiments of the present invention to illustrate the principles of the present invention, but are not intended to limit the scope of the present invention.

[0060] Example 1

[0061] A specific embodiment of the present invention discloses a multi-material topology optimization method based on a single-variable discrete interpolation scheme, such as... Figure 1 and Figure 2 As shown, it includes the following steps:

[0062] S1. Discretize the design domain using finite element methods and apply boundary conditions and external loads.

[0063] Specifically, finite element discretization is performed by meshing the design domain to obtain each finite element.

[0064] Specifically, boundary conditions define the physical behavior on the model boundaries; external loads refer to external forces or influences applied to the design domain.

[0065] S2. Based on the single-variable discrete combination interpolation method, a multi-material topology optimization model is constructed, and then the sensitivity of the multi-material topology optimization model is obtained.

[0066] In implementation, based on the univariate discrete combination interpolation method, a multi-material topology optimization model is constructed, including:

[0067] S21. Based on univariate characteristic functions, design variables are mapped to the topological density field of finite element elements. That is, based on univariate characteristic functions, a single design variable is mapped to multiple topological densities, and the topological density field is represented by different topological density functions.

[0068] In practical implementation, the univariate characteristic function is expressed as:

[0069]

[0070] In the formula, This indicates that when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... Topological density of phase candidate materials, This indicates that when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... The Heaviside mapping function for the candidate material, where η represents the penalty function. Θ Indicates the phase number of the candidate material.

[0071] Specifically, when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... Heaviside mapping function of phase candidate materials Represented as:

[0072]

[0073] In the formula, Indicates the first The threshold parameters for the phase candidate material are defined by β, where β represents the projection parameter, tanh() represents the hyperbolic tangent function, and H() represents the Heaviside function. The projection parameter β controls the steepness of the Heaviside function.

[0074] Specifically, the design variable χ(ξ) of the ξ-th finite element represents a value between 0 and 1.

[0075] The proposed method constructs multiple characteristic functions through a single optimization variable. These characteristic functions are associated with multiple discrete values ​​to distinguish different material phases. The design variable is transformed into a topological density field of multiple candidate materials using univariate characteristic functions.

[0076] S22. Based on the topological density field of the finite element, the filtered topological density field of the finite element is obtained, and then the normalized topological density field of the finite element is obtained.

[0077] Specifically, the topological density field of the filtered finite element is obtained in the following way:

[0078] First, the Helmholtz partial differential equation is expressed as:

[0079]

[0080] In the formula, Let ∠ denote the filtering radius, and ▽ denote the Laplace differential operator. This indicates that after filtering, when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element is... Topological density of candidate phase materials.

[0081] Then, based on formula (3), the following relationship is obtained after finite element discretization:

[0082]

[0083] In the formula, This represents the unfiltered density vector after discretization. Represented as the discretized filtered density vector, K ft T represents the filter matrix. ft This represents the transformation matrix.

[0084] Specifically, the filter matrix K ft and transformation matrix T ft Represented as:

[0085]

[0086] in,

[0087]

[0088]

[0089] In the formula, Ω represents the integration domain of the finite element element, and Sym. indicates that the matrix is ​​symmetric along the diagonal. These represent the x, y, and z coordinates of the first element node of the finite element element, respectively. These represent the x, y, and z coordinates of the second element node of the finite element element, respectively. These represent the x, y, and z coordinates of the third element node of the finite element element, respectively. These represent the x, y, and z coordinates of the fourth element node of the finite element element, respectively.

[0090] It should be noted that a coordinate system, such as a Cartesian coordinate system, is established based on the characteristics of the design domain and the requirements of the optimization problem. The origin of the coordinate system can be a specific point in the design domain, such as a symmetrical point or a mass point, and the direction of the coordinate axes can be determined according to the symmetry or other geometric features of the design domain.

[0091] Specifically, the filtered finite element topology density is mapped to the interval (0,1) using the Heaviside projection function, resulting in the normalized topology density field of the finite element, expressed as:

[0092]

[0093] In the formula, This indicates that after filtering, when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element is... Normalized topological density of phase candidate materials.

[0094] S23. Based on the normalized topological density field of the finite element, the elastic modulus of the candidate material is discretized and interpolated to obtain the elastic modulus of the finite element, and then the elastic tensor of the finite element is obtained.

[0095] Specifically, the elastic modulus of the candidate materials is discretized using the DMO method. The elastic modulus and elastic tensor of each finite element are weighted sums of the candidate materials, expressed as:

[0096]

[0097] In the formula, E[χ(ξ)] and D[χ(ξ)] represent the elastic modulus and elastic tensor of the ξ-th finite element when the design variable is χ(ξ), respectively; Let these represent the elastic modulus and elastic tensor of the empty material, respectively. They represent the first The elastic modulus and elastic tensor of the candidate materials, It indicates the ξ-th finite element. The weighting coefficients of candidate materials.

[0098] More specifically, the first elastic modulus of phase candidate materials Provided according to the specific application scenario.

[0099] More specifically, the first elastic tensor of candidate materials Represented as: According to the first elastic modulus of phase candidate materials The calculation of Poisson's ratio is as follows:

[0100] If it is a two-dimensional topology optimization, then

[0101]

[0102] If it is a three-dimensional topology optimization, then

[0103]

[0104] In the formula, It is the set number Poisson's ratio of the candidate material.

[0105] More specifically, to prevent singularities, the elastic modulus of hollow materials... It is usually set to a very small value, for example Then it can be based on The above formula is used for calculation.

[0106] More specifically, the ξ-th finite element element... Weighting coefficients of candidate materials Represented as:

[0107]

[0108] In the formula, Φ k [χ(ξ)] represents the design variables of the ξ-th finite element when χ(ξ) is the first finite element. Normalized topological density of k-phase candidate materials.

[0109] More specifically, the elastic modulus of the candidate material is an inherent material property of the candidate material, given according to the actual application scenario.

[0110] S24. Construct a multi-material topology optimization model based on the elastic tensor of finite element elements.

[0111] Specifically, based on the elastic tensor of the finite element element, the objective function and constraint function are established using finite element analysis technology to construct a multi-material topology optimization model.

[0112] Optionally, the multi-material topology optimization model is an optimization model with the objective function of minimizing structural flexibility and the constraint function of candidate material volume fraction, expressed as:

[0113]

[0114] In the formula, n represents the total number of finite element elements, C[u,χ(ξ)] represents the structural compliance of the ξ-th finite element under the structural displacement field u with the design variable value χ(ξ), ε(u) represents the linear strain tensor under the structural displacement field u, T represents the transpose, and a(u,δu) represents the structural compliance of the ξ-th finite element under the structural displacement field u and the Sobolev space H. 1 The bilinear properties of the virtual displacement field δu within (Ω), where l(δu) represents the bilinear properties in the Sobolev space H 1 The linear load function under the virtual displacement field δu within (Ω), where δu is the Sobolev space H 1 The virtual displacement field within (Ω), |Ω| represents the volume of the design domain Ω, Γ DThis represents the displacement boundary, i.e., the Dirichlet boundary condition, where g represents the condition at the displacement boundary Γ. D displacement vector on, Indicates the first Volume constraint function for candidate phase materials, Indicates the first The maximum volume fraction of the candidate material, υ0 represents the volume of each finite element.

[0115] Specifically, in the structural displacement field u and the Sobolev space H 1 The bilinear property a(u,δu) under the virtual displacement field δu within (Ω) is expressed as:

[0116]

[0117] In the formula, ε(δu) represents the linear strain tensor under the virtual displacement field δu.

[0118] Specifically, in the Sobolev space H 1 The linear load function l(δu) under the virtual displacement field δu within (Ω) is expressed as:

[0119]

[0120] In the formula, f represents the external force acting on the structure, and h represents the force at the boundary Γ. N The traction force distributed on the top, Γ N This represents the force boundary, i.e., the Neumann boundary condition. Here, the external forces acting on the structure and the forces at the boundary Γ represent the forces acting on the structure. N The traction force distributed on the top is set according to the actual application scenario.

[0121] Correspondingly, in the optimization model with minimizing structural flexibility as the objective function and candidate material volume fraction as the constraint function, the sensitivities of the objective function and the constraint function are expressed as follows:

[0122]

[0123] in,

[0124]

[0125] In the formula, K ft T represents the filter matrix. ft The transformation matrix is ​​represented by sech(), which represents the hyperbolic secant function. Represents the first of the matrix All elements of the column, This indicates that when the design variable value of the ξ-th finite element is χ(ξ), the ξ-th finite element... Topological density of candidate phase materials.

[0126] Optionally, the multi-material topology optimization model is an optimization model with the objective function of minimizing structural flexibility and the constraint function of candidate material mass density, expressed as:

[0127]

[0128] In the formula, G m M represents the mass density constraint function. f Indicates the maximum mass density. Indicates the first The mass density of the phase candidate materials. The mass density of each phase candidate material is given according to the actual application scenario.

[0129] Correspondingly, for the optimization model with the objective function of minimizing structural flexibility and the candidate material mass density as the constraint function, the sensitivities of the objective function and the constraint function are expressed as follows:

[0130]

[0131] It should be noted that the sensitivity of the multi-material topology optimization model is derived as follows:

[0132] First, the derivative of the objective function with respect to the design variable χ(ξ) is derived using the chain rule:

[0133]

[0134] Wherein, the objective function has weight coefficients W ξ,j The sensitivity is expressed as:

[0135]

[0136] According to the rule of differentiation, we can obtain:

[0137]

[0138] Substituting equation (23) into equation (22), we get:

[0139]

[0140] The first derivative of the elastic tensor with respect to the weighting coefficients is obtained from equation (8). Substituting this into equation (24), we get:

[0141]

[0142] Based on equation (11), the weight coefficients are obtained. Sensitivity:

[0143]

[0144] Based on equation (6), we obtain about The partial derivatives are expressed as:

[0145]

[0146] Combining the Helmholtz partial differential filtering equations given in equations (3) and (4), the filtered topological density field is obtained. Topological density field before filtering The sensitivity is as follows:

[0147]

[0148] According to equations (1) and (2), the partial derivatives of the univariate characteristic function with respect to the design variables are expressed as:

[0149]

[0150] Substituting equations (25), (26), (27), (28), and (29) into equation (21), we obtain the sensitivity of the objective function to the design variables:

[0151]

[0152] in,

[0153]

[0154] In the formula, K ft T represents the filter matrix. ft The transformation matrix is ​​represented by sech(), which represents the hyperbolic secant function. This indicates that when the design variable value of the ξ-th finite element is χ(ξ), the ξ-th finite element... Topological density of candidate phase materials.

[0155] Then, the sensitivity of the constraint function to the design variables is derived.

[0156] According to the chain rule, the volume constraint function The sensitivity to design variables can be expressed as follows:

[0157]

[0158] The volume constraint can be obtained from the volume constraint function given by equation (12). right The partial derivatives are shown below:

[0159]

[0160] Substituting equations (27), (28), (29), and (33) into equation (32), we obtain the sensitivity of the volume constraint function to the design variables:

[0161]

[0162] Similarly, the mass constraint function G can be given. m Sensitivity to design variable χ(ξ):

[0163]

[0164] S3. Based on the multi-material topology optimization model and sensitivity, update the design variables to obtain the optimal topology density field.

[0165] During implementation, based on the multi-material topology optimization model and sensitivity, the gradient optimization algorithm MMA is used to update the design variables until the multi-material topology optimization model converges, thereby obtaining the optimal design variables and the optimal topological density field. These are then substituted into the post-processing function to draw the material distribution contour map.

[0166] In practice, when updating design variables using the gradient optimization algorithm MMA, the optimization parameters in the iteration are updated based on an adaptive parameter adjustment strategy; the optimization parameters include penalty parameters and projection parameters.

[0167] Specifically, the adaptive parameter adjustment strategy is as follows:

[0168] Before the gradient optimization algorithm MMA iteration, the optimization parameters are initialized, and the value of the counter a is initialized to 0;

[0169] During the MMA gradient optimization algorithm iteration, it is determined whether the projection parameter β′ in the current iteration is less than the set projection parameter threshold. If so, then the design variables in the current iteration are taken as the optimal design variables; otherwise:

[0170] Determine if the counter value is greater than 0 in the current iteration.

[0171] If the value of the counter 'a' in the current iteration is less than or equal to 0, then determine whether the change δ of the objective function in the current iteration compared to the objective function in the previous iteration is less than a predetermined value for the change. If yes, then set the counter value a to 1, and do not update the optimization parameters in the current iteration; if no, do not update the optimization parameters in the current iteration.

[0172] If the value of the counter in the current iteration is greater than 0, the value of the counter a is incremented by 1, and it is determined whether the value of the counter a is at the set number of iterations threshold. If it is, the optimization parameters in the current iteration are updated, and the value of the counter a is set to 1; otherwise, the optimization parameters in the current iteration are not updated.

[0173] More specifically, the projection parameters in the optimization parameters are updated using the following formula:

[0174] β=θ×β′ (36)

[0175] In the formula, β′ is the projection parameter in the current iteration, β is the updated projection parameter, and θ represents the first adjustment coefficient.

[0176] Specifically, the penalty parameter in the optimization parameters is updated using the following formula:

[0177]

[0178] In the formula, η′ is the penalty parameter in the current iteration, η is the penalty parameter after the update, and σ represents the second adjustment coefficient.

[0179] More specifically, the predetermined value of the change is set to 1e-4; the threshold value of the counter is set to 4.

[0180] More specifically, the first and second adjustment coefficients are set according to actual requirements. The projection parameter thresholds are set according to actual needs.

[0181] More specifically, in the optimization parameters, the penalty parameter η is initialized to 1; the projection parameter β is initialized to a small positive value; for example, the projection parameter β = 0.3.

[0182] Specifically, based on the obtained optimal design variables, the optimal topological density field of each phase candidate material is further obtained through formula (11).

[0183] It is understandable that the optimization method in this embodiment can adaptively adjust the penalty parameter and projection parameter, thereby making the algorithm more robust and adaptable during the iteration process, avoiding local optima and numerical instability; it can achieve a reasonable distribution of materials in the structure while ensuring material mass or volume fraction constraints; in the optimization design, low-stiffness materials are mainly distributed in non-critical areas, medium-stiffness materials are distributed on the load transfer path, and high-stiffness materials are distributed in the load application area and support area, with a clear load transfer path; in some structures, the use of intermediate materials can be avoided.

[0184] Compared with existing technologies, this invention provides a multi-material topology optimization method based on a single-variable discrete interpolation scheme. By discretizing the design domain using finite element methods and applying boundary conditions and external loads, a multi-material topology optimization model is constructed based on a single-variable discrete combination interpolation method. The sensitivity of the multi-material topology optimization model is then obtained. Based on the multi-material topology optimization model and the sensitivity, the design variables are updated to obtain the optimal topological density field. This solves the problems of numerous design variables, low optimization efficiency, unreasonable design results, and difficulty in manufacturing the designed structure. By constructing a single-variable characteristic function and using a single variable to select multiple discrete physical quantities, the composite function of the reduced-order multi-characteristic function of the single variable is combined with the corresponding material properties for material interpolation. This eliminates the material envelope phenomenon present in existing step-type interpolation models and overcomes the disadvantage of a limited design domain. The number of design variables is independent of the number of candidate material phases, significantly reducing the number of design variables. While achieving similar or even lower objective function values, fewer iterations are needed, improving optimization efficiency. It can achieve a reasonable allocation of materials in the structure while ensuring material mass or volume fraction constraints.

[0185] Example 2

[0186] To verify the effectiveness and correctness of the method proposed in Example 1, specific examples are provided for verification.

[0187] Scenario 1: Using a volume constraint approach, material allocation is performed by setting the volume fraction of the design domain occupied by candidate materials to minimize structural flexibility, thus verifying the effectiveness of the proposed two-dimensional topology optimization method.

[0188] like Figure 3 As shown, taking a two-dimensional MBB beam as an example, the parameters are set to length unit L = 60 and external force F = 1; half of the structure is selected as the design domain, and a 120×60 quadrilateral grid is used for structural discretization; the elastic parameters of the candidate materials and the given volume fraction constraints are shown in Table 1; the first adjustment coefficient of the projection parameter β is θ = 1.24, and the second adjustment coefficient of the penalty parameter η is σ = 1.065.

[0189] Table 1 Elastic parameters of candidate materials and given volume fraction constraints

[0190]

[0191] like Figure 4 The figure shows the objective function value and the iterative convergence process of the volume fraction in the MBB beam structure. As can be seen from the figure, the objective function value first rises and then falls rapidly in the early stage of iteration, and then tends to stabilize after about 50 iterations; the volume fraction curves show that the volume fractions of MAT1 and MAT2 also gradually converge during the iteration process, and finally approach the set volume constraint of 0.25 respectively.

[0192] Figures 5(a) and 5(b) show the comparison between the two-dimensional MBB beam topology optimization results obtained by the method proposed in Example 1 and the DMO method under volume constraints, for the two cases θ = 1.24, σ = 1.065 and θ = 1.1445, σ = 1.0247 respectively. Specifically, when θ = 1.24 and σ = 1.065, the DMO method has 89 iterations, an objective function of 23.9965, a volume fraction of 0.2499 for material 1, and a volume fraction of 0.2498 for material 2; when θ = 1.24 and σ = 1.065, the method in Example 1 has 88 iterations, an objective function of 23.7622, ​​a volume fraction of 0.2497 for material 1, and a volume fraction of 0.2497 for material 2. When θ = 1.1445 and σ = 1.0247, the number of iterations in the DMO method is 189, the objective function is 24.0376, the volume fraction of material 1 is 0.2500, and the volume fraction of material 2 is 0.2500. When θ = 1.1445 and σ = 1.0247, the number of iterations in the method of Example 1 is 191, the objective function is 23.5106, the volume fraction of material 1 is 0.2500, and the volume fraction of material 2 is 0.2500.

[0193] As can be seen from Figures 5(a) and 5(b), under the premise of satisfying the volume fraction constraint, the method proposed in Example 1 achieves a lower objective function compliance value with a similar number of iterations as the DMO method, indicating that its optimization effect is better than that of the DMO method. It should be noted that the material distribution obtained by the traditional DMO method is relatively mixed, increasing the manufacturing difficulty. Conversely, the material distribution obtained by the method proposed in Example 1 is concentrated, with hard materials allocated to areas with higher stress and soft materials allocated to the middle areas, resulting in an optimized structure that is more suitable for force transmission and load bearing. Furthermore, for the first and second adjustment coefficients, the smaller their values, the better the objective function value will be, but a higher number of iterations is also required.

[0194] Scenario 2: Using a mass constraint approach, material allocation is performed by setting the mass of candidate materials to minimize structural flexibility, thus verifying the effectiveness of the proposed method in two-dimensional topology optimization.

[0195] like Figure 6 As shown, taking a cantilever beam structure as an example, the length unit L = 80, the design domain is discretized by a 120×80 quadrilateral element with four nodes, and the middle of the right segment bears a unit load F.

[0196] Consider the multi-material topology optimization design of the cantilever beam under a single mass constraint. There are five materials, and the material parameters are shown in Table 2. The mass constraint is set to 0.3.

[0197] Table 2 Material Parameters

[0198]

[0199] like Figure 7 As shown, the objective function and volume fraction change with the number of iterations for the five material quantities. It can be seen that the objective function converges after 60 iterations. MAT1 and MAT2 also remain unchanged after 60 iterations, while MAT3, MAT4, and MAT5 converge after 50 iterations. It is worth noting that the volume fraction of MAT4 gradually decreases, eventually converging to 0.

[0200] like Figure 8(a) and 8(b) The diagram shows a comparison between the method proposed in Example 1 and the DMO method when using five materials (MAT1 / 2 / 3 / 4 / 5) and three materials (MAT2 / 3 / 5). Specifically, when using five materials (MAT1 / 2 / 3 / 4 / 5), the DMO method requires 82 iterations, has an objective function of 39.2550, and the volume fractions of material 1, 2, 3, 4, and 5 are 0.1266, 0.3500, 0.1250, 0.0032, and 0.0031, respectively. In Example 1, the method requires 86 iterations, has an objective function of 36.8620, and the volume fractions of material 1, 2, 3, 4, and 5 are 0.5576, 0.1869, 0.0544, 0.0000, and 0.0037, respectively. When using three materials (MAT2 / 3 / 5), the DMO method has 87 iterations, an objective function of 39.3614, a volume fraction of 0.2487 for material 1, a volume fraction of 0.2614 for material 2, and a volume fraction of 0.0056 for material 3. The method proposed in Example 1 has 86 iterations, an objective function of 38.2253, a volume fraction of 0.5324 for material 1, a volume fraction of 0.0073 for material 2, and a volume fraction of 0.0289 for material 3.

[0201] As can be seen, for the topology optimization results under the five materials, the optimal objective function value obtained by the DMO method is 39.2550, while the optimal objective function value obtained by the proposed method is 36.8620, which is significantly better than the DMO method. In the final topology optimization configuration, material MAT4 is completely degraded, and material MAT5 with the highest stiffness is selected to fill the load application area, while material MAT3 is filled near the support end. It should be noted that the material distribution obtained by the DMO method is relatively scattered, while the material distribution obtained by the proposed method is concentrated. For the topology optimization results under the three materials, the minimum compliance value of the DMO method is 39.3614, while the minimum compliance value of the method proposed in Example 1 is 38.2253, which also shows better performance. In terms of material distribution, the method in Example 1 fills the load location and the support end location with MAT5, and the material distribution is more reasonable.

[0202] Scenario 3: Using a volume constraint approach, material allocation is performed by setting the volume fraction of the design domain occupied by candidate materials to minimize structural flexibility, thus verifying the effectiveness of the proposed method in three-dimensional topology optimization.

[0203] Taking a bridge design problem as an example, the structural tetrahedral mesh was generated using the open-source 3D finite element mesh generator GMSH. The candidate materials were aluminum, iron, and titanium, with volume fractions of 0.2, 0.15, and 0.1, respectively. The material parameters are shown in Table 3.

[0204] Table 3 Elastic parameters of candidate materials for 3D bridge structures

[0205]

[0206] like Figure 9 As shown, the loads and boundary conditions for the bridge design problem are presented. The overall dimensions of the bridge structure are 400mm × 80mm × 80mm. In this example, the bridge structure is divided into a design domain Ω. D Non-design domain Ω ND Multi-material topology optimization is performed by assigning design variables within the design domain, while no optimization variables are assigned outside the design domain, where material distribution and structural configuration remain unchanged. A uniformly distributed force F = 16667 Pa is applied along ξ. z The negative direction of the axis acts on the non-design area Ω of the bridge ND Above. Fixed displacement areas Ω at both ends below the bridge. B The boundary conditions are defined. The entire domain has 63,463 grid cells and 14,267 nodes. Based on the proposed multi-material topology optimization method, candidate materials are rationally allocated at different spatial locations of the bridge to minimize the flexibility of the bridge structure while satisfying volume constraints. Due to the symmetry of the structure, one-quarter of the structure is used during the optimization process. In post-processing, the topology and material distribution of the complete structure are obtained through mirroring.

[0207] The objective function value and iterative convergence of the volume fractions of the three materials in the multi-material topology optimization process of the bridge structure are as follows: Figure 10 As shown, the objective function value initially rises rapidly during the initial iterations, then decreases rapidly after the 40th iteration, and gradually stabilizes. Furthermore, the volume fraction of aluminum decreases rapidly in the early stages of iteration and then stabilizes. The other two materials show a trend of first decreasing and then increasing. Around the 50th iteration, the volume fractions of the three materials essentially converge. This indicates that during the optimization process, the material distribution gradually reaches equilibrium, minimizing the objective function and verifying the effectiveness and stability of the multi-material topology optimization method.

[0208] like Figure 11 As shown in the figure, the material distribution results of the three-dimensional bridge structure after multi-material topology optimization are presented. Through the topology optimization results displayed from different perspectives in the figure, it can be seen that the high-stiffness metallic material titanium is concentrated in the critical stress area and the support area. The medium-stiffness metallic material iron is mainly distributed along the load transfer path, while the low-stiffness metallic material aluminum is mainly distributed in non-critical areas, playing a connecting role. The obtained bridge structure topology and material distribution not only improve the structural stiffness but also enhance the load transfer path, and the resulting material distribution boundaries are clear and smooth. Therefore, it can be concluded that the multi-material topology optimization method proposed in this paper can effectively achieve reasonable material allocation in three-dimensional structures to minimize structural flexibility, verifying the feasibility of this method in the multi-material topology optimization design of complex three-dimensional structures.

[0209] Figure 12 A comparison of the results of 3D bridge structure topology optimization using the traditional DMO method and the method proposed in Example 1 is presented. The DMO method has 89 iterations, an objective function of 1.7675e-4, and the volume fractions of material 1, 2, and 3 are 0.1998, 0.1498, and 0.1000, respectively. The method proposed in Example 1 has 79 iterations, an objective function of 1.4198, and the volume fractions of material 1, 2, and 3 are 0.2032, 0.1494, and 0.0999, respectively.

[0210] It can be seen that the results obtained by the method proposed in Example 1 are largely consistent with the topology obtained by the DMO method. Regarding material distribution, both methods distribute high-stiffness titanium material in the loaded region. The difference lies in the distribution method: Example 1 distributes iron in the bridge crown region, while the DMO method distributes aluminum. Furthermore, Example 1 distributes high-stiffness titanium in the boundary fixed support region, while the DMO method also distributes aluminum. From the perspective of load transfer path, the material distribution obtained by the method proposed in Example 1 is more favorable. It should be noted that the method proposed in Example 1 achieves a lower objective function value with fewer iterations. This demonstrates the advantage of the method proposed in Example 1.

[0211] Those skilled in the art will understand that all or part of the processes of the methods described in the above embodiments can be implemented by a computer program instructing related hardware, and the program can be stored in a computer-readable storage medium. The computer-readable storage medium may be a disk, optical disk, read-only memory, or random access memory, etc.

[0212] The above description is only a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in the present invention should be included within the scope of protection of the present invention.

Claims

1. A multi-material topology optimization method based on a univariate discrete interpolation scheme, characterized in that, Includes the following steps: The design domain is discretized using finite element methods, and boundary conditions and external loads are applied. Based on the univariate discrete combination interpolation method, a multi-material topology optimization model is constructed, and then the sensitivity of the multi-material topology optimization model is obtained; Based on the multi-material topology optimization model and sensitivity, the design variables are updated to obtain the optimal topology density field; The construction of the multi-material topology optimization model based on the univariate discrete combination interpolation method includes: Based on univariate characteristic functions, design variables are mapped to the topological density field of finite element elements; Based on the topological density field of the finite element, the filtered topological density field of the finite element is obtained, and then the normalized topological density field of the finite element is obtained. Based on the normalized topological density field of the finite element, the elastic modulus of the candidate material is discretized and interpolated to obtain the elastic modulus of the finite element, and then the elastic tensor of the finite element is obtained. A multi-material topology optimization model is constructed based on the elastic tensor of finite element elements.

2. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 1, characterized in that, The univariate feature function is expressed as: In the formula, This indicates that when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... Topological density of phase candidate materials, This indicates that when the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element... The Heaviside mapping function for the candidate material is given by η, where η represents the penalty function and Θ represents the number of phases of the candidate material.

3. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 2, characterized in that, When the design variable of the ξ-th finite element is χ(ξ), the ξ-th finite element is χ(ξ). Heaviside mapping function of phase candidate materials Represented as: In the formula, Indicates the first The threshold parameter for the candidate material, β represents the projection parameter, and tanh() represents the hyperbolic tangent function.

4. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 3, characterized in that, The elastic modulus and elastic tensor of the finite element are expressed as follows: In the formula, E[χ(ξ)] and D[χ(ξ)] represent the elastic modulus and elastic tensor of the ξ-th finite element when the design variable is χ(ξ), respectively; Let these represent the elastic modulus and elastic tensor of the empty material, respectively. They represent the first The elastic modulus and elastic tensor of the candidate materials, It indicates the ξ-th finite element. The weighting coefficients of candidate materials.

5. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 4, characterized in that, In the ξth finite element element Weighting coefficients of candidate materials Represented as: In the formula, Φ k [χ(ξ)] represents the design variables of the ξ-th finite element when χ(ξ) is the first finite element. Normalized topological density of k-phase candidate materials.

6. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 1, characterized in that, Based on the multi-material topology optimization model and sensitivity, the gradient optimization algorithm MMA is used to update the design variables until the multi-material topology optimization model converges, thereby obtaining the optimal design variables and the optimal topology density field.

7. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 6, characterized in that, When updating design variables using the gradient optimization algorithm MMA, the optimization parameters in the iteration are updated based on an adaptive parameter adjustment strategy; the optimization parameters include penalty parameters and projection parameters.

8. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 7, characterized in that, The adaptive parameter adjustment strategy includes: Before the gradient optimization algorithm MMA iteration, the optimization parameters are initialized, and the counter value is initialized to 0; During the MMA gradient optimization algorithm iteration, it is determined whether the projected parameters in the current iteration are less than the set projection parameter threshold. If so, the design variable in the current iteration is taken as the optimal design variable; otherwise: Determine if the counter value is greater than 0 in the current iteration. If yes, then determine whether the change in the objective function in the current iteration is less than the predetermined value of the change. If yes, then set the counter value to 1 and do not update the optimization parameters in the current iteration; otherwise, do not update the optimization parameters in the current iteration. If not, the counter value is incremented by 1, and it is determined whether the current counter value is the set number of times threshold. If yes, the optimization parameters in the current iteration are updated, and the counter value is set to 1; otherwise, the optimization parameters in the current iteration are not updated.

9. The multi-material topology optimization method based on a single-variable discrete interpolation scheme according to claim 8, characterized in that, The projection parameters in the optimization parameters are updated using the following formula: β=θ×β′ In the formula, β′ is the projection parameter in the current iteration, β is the updated projection parameter, and θ represents the first adjustment coefficient; The penalty parameter in the optimization parameters is updated using the following formula: In the formula, η′ is the penalty parameter in the current iteration, η is the penalty parameter after the update, and σ represents the second adjustment coefficient.

Images