A Thermomechanical Coupling Topology Optimization Method for Thin-walled Non-uniform Filling Structures

Abstract

The present invention provides a method for thermomechanical coupling topology optimization of thin-walled non-uniform filling structures, step 1: establish a finite element model according to design requirements, and determine the non-design domain and the design domain; step 2: initialize the design variables, and then obtain the local volume fraction aggregation field through one-step filtering; step 3: describe the distribution of materials in the design domain; step 4: perform thermomechanical coupling finite element calculations to solve the structural displacement field U; step 5: perform corresponding sensitivity analysis; step 6: submit the sensitivity to the moving asymptote algorithm solver to update the overall topological field and lattice density field; step 7: check whether the iteration stop condition is met. The present invention realizes the lightweight design of the power battery module bracket structure for thermoelastic problems; the power battery module bracket with thin-walled non-uniform filling structure characteristics improves the structural stiffness and exhibits good heat dissipation performance.

Description

Technical Field

[0001] The invention relates to the technical field of ultra-lightweight structure design methods, and in particular to a thermal-mechanical coupling topology optimization method for a thin-walled non-uniformly filled structure. Background Art

[0002] The thin-walled non-uniform filling structure is a lightweight composite structure with excellent performance, which is mainly composed of an outer thin-walled layer and an internal non-uniform lattice material filling. Theoretical and research results show that the thin-walled filling structure has the advantages of high specific stiffness, large tensile strength, good anti-buckling performance and energy absorption. On the other hand, the non-uniform filling structure itself contains many small structures and voids, and also has great potential in heat conduction. Through the coordinated design of the overall structure and the filling configuration, the filling area structure can be optimized simultaneously, which can maximize the mechanical properties and lightweight level of the structure.

[0003] As a structural optimization method, topology optimization can maximize the potential of design freedom and often obtain novel structures in the conceptual design stage. Thermomechanical coupling structural topology optimization design can comprehensively consider the topology optimization design of the internal heat transfer structure of the structure, that is, by reasonably arranging high thermal conductivity materials inside the structure, designing the optimal heat transfer path to improve the heat dissipation capacity of the structure. At present, the structures of thermomechanical coupling topology optimization are often solid filling or uniform filling configurations, which cannot adapt to complex mechanical-thermal stress distribution and have limited lightweight potential. Thin-walled non-uniform filling structures can be obtained by using local volume constraints. Further thermomechanical coupling topology optimization can fully explore the lightweight potential of the structure while ensuring mechanical properties and heat dissipation performance. Summary of the invention

[0004] Based on local volume constraints, unidirectional steady-state thermal-mechanical coupling finite element analysis, two-step filtering and Heaviside projection method, and multi-scale material interpolation model of thin-walled filling structures, the present invention constructs a thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures, which is expected to provide technical support for the high performance and lightweight design of carrier structures.

[0005] Specifically, the present invention proposes a thermal-mechanical coupling topology optimization method for a thin-walled non-uniformly filled structure, and the thermal-mechanical coupling topology optimization method for a thin-walled non-uniformly filled structure comprises the following steps:

[0006] Step 1: Establish a finite element model according to the design requirements, determine the non-design domain and design domain, select the overall topology field μ and the lattice density field v as design variables, establish an optimization model with structural flexibility as the target and mass fraction as the constraint, and initialize the model including finite element discretization, definition of mechanical field and temperature field finite element calculation related parameters, PDE filter matrix preprocessing and iteration related parameter initialization;

[0007] Step 2: Initialize the design variable μ and obtain the density field through a two-step filtering algorithm and the normalized spatial gradient norm τ; initialize the design variable ν, perform one-step filtering and projection to obtain the non-uniform filling density field ψ, and then perform one-step filtering to obtain the local volume fraction aggregation field

[0008] Step 3: By converting the density field Spatial gradient norm τ and local volume fraction aggregation field Substitute the material interpolation function to calculate the elastic modulus E, thermal conductivity k, and thermal stress coefficient f β and mass density ρ to describe the distribution of material in the design domain;

[0009] Step 4: Assemble the overall thermal conductivity matrix K th And the overall stiffness matrix K, perform thermal-mechanical coupling finite element calculation to solve the structural displacement field U;

[0010] Step 5: Calculate the structural flexibility, global mass fraction, local volume constraint and node temperature constraint values, and perform corresponding sensitivity analysis;

[0011] Step 6: Determine the objective function and constraints according to the mathematical model of the optimization problem, and submit the sensitivity to the moving asymptote algorithm solver to update the overall topological field μ and the lattice density field v;

[0012] Step 7: Check whether the iteration stop condition is met: the maximum change of the design variable is less than 0.01. If so, output the result; otherwise, return to step 2 and repeat the above steps.

[0013] Furthermore, in step 1, the optimization model can be expressed as:

[0014]

[0015] ikB th T=Q

[0016] KU=F m +F th

[0017]

[0018] lvc(v)≤0

[0019] gmc(μ,v)=M(μ,v) / M * -1≤0

[0020] 0≤μ e ≤1,0≤v e ≤1

[0021] In the formula, st represents the constraint condition, T represents the transpose of the matrix, c represents the flexibility, and μ represents the relative density of the unit. e The vector composed of v is the relative density of the unit v e The vector composed of K is the overall stiffness matrix of the structure, K th is the overall heat transfer matrix of the structure, U and T are the node displacement vector and node temperature vector respectively, Q is the node heat load vector, F m and F th are the node mechanical load vector and the equivalent node thermal load vector, t i is the temperature value of the ith node, is the temperature constraint value of the ith node, lvc is the local volume constraint of the non-uniform filling area, and gmc is the global mass constraint.

[0022] Furthermore, in step 2, the filtering expression is:

[0023]

[0024] In the formula, is the Laplace operator, τ is the cell density before filtering, and the cell density after filtering can be obtained by applying the filtering boundary conditions. R f is the filter radius.

[0025] Furthermore, in step 2, the cell density after filtering of the e-th cell is The relative density of the unit obtained after Heaviside projection for:

[0026]

[0027] Where β is the sharpness of the projection, η is the threshold of the projection, Represents the relative density of the e-th unit after filtering The relative density of the unit obtained after projection;

[0028] Through the two-step filtering and projection operation on μ, the filling density field can be obtained and the thin-wall density field By one-step filtering and projection operation on v, the non-uniform filling density field ψ can be obtained. Applying a smoothing filter on the filling field ψ can obtain the local volume fraction field of the structure:

[0029] Furthermore, in step 3, the material interpolation function is:

[0030]

[0031] Where ρ, E, κ are the overall mass density, elastic modulus and thermal conductivity, ρ 0 、E 0 and κ 0 They represent the mass density, elastic modulus and thermal conductivity of the shell region material, respectively, ρ , E and λ κ They represent the ratio of the material properties of the filling region to the mass density, elastic modulus and thermal conductivity of the shell region material, respectively. In order to avoid numerical singularity in the finite element calculation of the empty region material, the minimum quantity ρ is introduced. min 、E min and κ min .

[0032] Similarly, the thermal stress coefficient f described in step 3 β The interpolation function is as follows:

[0033]

[0034] In the formula, E 0 represents the elastic modulus of the shell region material, λ E It represents the ratio of the material properties of the filling area to the elastic modulus of the shell area material. In order to avoid numerical singularity in the finite element calculation of the empty phase area material, the minimum quantity E is introduced. min .

[0035] Furthermore, in step 4, the overall thermal conductivity matrix K th And the global stiffness matrix K is assembled as follows:

[0036]

[0037] In the formula, ∑ represents the overall matrix assembly operator, k 0 and is the unit stiffness matrix and the unit heat conduction matrix.

[0038] Furthermore, in step 5, the sensitivity of the flexibility c to the design parameters can be obtained by the following formula:

[0039]

[0040] Where E is the elastic modulus, k is the thermal conductivity, and f β are thermal stress coefficients, both of which are related to the design variables μ and v.

[0041] Furthermore, in step 5, the sensitivity of the global mass constraint gmc with respect to the design parameters can be obtained by the following formula:

[0042]

[0043] According to the chain rule, the derivative of the improved mass density ρ′ with respect to the design variable μ is:

[0044]

[0045] Since the improved mass density ρ′ is independent of the design variable v,

[0046]

[0047] Furthermore, in step 5, the node temperature constraint sensitivity and the local volume constraint sensitivity can be obtained by the following formula:

[0048]

[0049] The beneficial effects achieved by the present invention are:

[0050] 1. The present invention adopts a two-step filtering and projection method to identify the characteristics of thin-walled filling structures. Through local volume constraints, the design variables are driven to generate non-uniform filling structures with micro-grids. Combining the thin-wall density field and the non-uniform filling density field, an interpolation model of thin-walled filling structure materials with non-uniform filling characteristics is further constructed.

[0051] 2. The present invention takes structural flexibility as the optimization target, adopts unidirectional thermomechanical coupling finite element analysis, and takes mass fraction as the constraint to construct a thermomechanical coupling topology optimization model for thin-walled non-uniform filling structures; the adjoint sensitivity analysis method is used to derive the objective function, node temperature constraints, local volume constraints and global mass fraction constraints, thereby achieving robust thermomechanical coupling topology optimization of thin-walled filling structures.

[0052] 3. The present invention utilizes the proposed thin-walled non-uniform filling structure thermo-mechanical coupling topology optimization method to achieve a lightweight design of a power battery module bracket structure for thermal elastic problems; the power battery module bracket with thin-walled non-uniform filling structure characteristics improves the structural stiffness while exhibiting good heat dissipation performance. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] Figure 1 A flow chart of a thermal-mechanical coupling topology optimization method for a thin-walled non-uniform filling structure constructed for the present invention;

[0054] Figure 2 A schematic diagram of a two-dimensional power battery module support structure design domain, loads and boundary conditions provided in an embodiment of the present invention;

[0055] Figure 3 This is a list of initial design parameters for a thin-walled non-uniform filling structure power battery module bracket in an embodiment of the present invention;

[0056] Figure 4A schematic diagram of a thin-walled non-uniform filling density field described by the present invention through two-step PDE filtering and projection operations;

[0057] Figure 5 Schematic diagram of an optimized thin-wall non-uniform filling structure in an embodiment of the present invention;

[0058] Figure 6 Schematic diagram of temperature field distribution of the optimized thin-walled non-uniform filling structure in an embodiment of the present invention.

[0059] Figure 7 Schematic diagram of the optimized power battery module support structure with a thin-walled non-uniform filling grid structure in an embodiment of the present invention. DETAILED DESCRIPTION

[0060] The technical solution of the present invention is described in more detail below in conjunction with the accompanying drawings. The present invention includes but is not limited to the following embodiments.

[0061] This example presents a thermal-mechanical coupling topology optimization method for thin-walled non-uniformly filled structures. The design process is as follows: Figure 1 As shown, the specific steps include:

[0062] Step 1: Establish a finite element model according to the design requirements, determine the non-design domain and design domain, select the overall topology field μ and the lattice density field v as design variables, establish an optimization model with structural flexibility as the target and mass fraction as the constraint, and initialize the model including finite element discretization, definition of mechanical field and temperature field finite element calculation related parameters, PDE filter matrix preprocessing and iteration related parameter initialization;

[0063] In this implementation, the established optimization model can be expressed as:

[0064]

[0065] ikB th T=Q

[0066] KU=F m +F th

[0067]

[0068] lvc(v)≤0

[0069] gmc(μ,v)=M(μ,v) / M * -1≤0

[0070] 0≤μ e ≤1,0≤v e ≤1

[0071] In the formula, st represents the constraint condition, T represents the transpose of the matrix, c represents the flexibility, and μ represents the relative density of the unit. e The vector composed of v is the relative density of the unit v e The vector composed of K is the overall stiffness matrix of the structure, K th is the overall heat transfer matrix of the structure, U and T are the node displacement vector and node temperature vector respectively, Q is the node heat load vector, F m and F th are the node mechanical load vector and the equivalent node thermal load vector, t i is the temperature value of the ith node, is the temperature constraint value of the ith node, lvc is the local volume constraint of the non-uniform filling area, and gmc is the global mass constraint.

[0072] Step 2: Initialize the design variable μ and obtain the density field through a two-step filtering algorithm and the normalized spatial gradient norm τ; initialize the design variable ν, perform one-step filtering and projection to obtain the non-uniform filling density field ψ, and then perform one-step filtering to obtain the local volume fraction aggregation field

[0073] First, filtering is used to avoid the checkerboard phenomenon and grid dependence, and it helps to obtain spatial gradient information. The solution format is as follows:

[0074]

[0075] In the formula, is the Laplace operator, τ is the cell density before filtering, and the cell density after filtering can be obtained by applying the filtering boundary conditions. R f is the filter radius.

[0076] Based on filtering, combined with Heaviside projection, a clear 0-1 distribution topological configuration can be obtained. The cell density after filtering of the e-th cell is The relative density of the unit obtained after Heaviside projection for:

[0077]

[0078] Where β is the sharpness of the projection, η is the threshold of the projection, Represents the relative density of the e-th unit after filtering The relative density of the unit obtained after projection;

[0079] Through two-step filtering and projection operations, the filling density field can be obtained Thin wall density field Inhomogeneous filling density field ψ and local volume fraction field of the structure

[0080] Specifically, the design variable μ is filtered for the first time with R1 as the radius to obtain the unit density Again Projection is performed to obtain the filling density field Then take R2 as the radius Perform a second filter to obtain the density field beg The spatial gradient norm field ( is the gradient operator), and multiplied by the normalization factor The normalized spatial gradient norm field can be obtained Then perform a second projection to obtain

[0081] In addition, the design variable v is filtered for the first time with R3 as the radius to obtain the unit density Again Projection to obtain a non-uniform filling density field The local volume fraction field of the structure is obtained by filtering ψ for the second time with R4 as the radius

[0082] Step 3: By converting the density field Spatial gradient norm τ and local volume fraction aggregation field Substitute the material interpolation function to calculate the elastic modulus E, thermal conductivity κ, and thermal stress coefficient f β and mass density ρ to describe the distribution of material in the design domain;

[0083] Specifically, in order to obtain the material model describing the distribution of the non-uniformly filled sandwich structure in the design domain, ρ, E, and κ are used to represent the overall mass density, elastic modulus, and thermal conductivity, respectively. 0 、E 0 and κ 0 They represent the mass density, elastic modulus and thermal conductivity of the shell region material, respectively, ρ , E and λ κ They represent the ratio of the material properties of the filling region to the mass density, elastic modulus and thermal conductivity of the shell region material. In addition, in order to avoid numerical singularities in the finite element calculation of the empty region material, the minimum quantity ρ is introduced. min 、E min and k min . The material interpolation function can be established as follows:

[0084]

[0085] Similarly, the interpolation function of the thermal stress coefficient can be obtained as follows:

[0086]

[0087] Step 4: Assemble the overall thermal conductivity matrix K th And the overall stiffness matrix K, perform thermal-mechanical coupling finite element calculation to solve the structural displacement field U;

[0088] Specifically, the overall thermal conductivity matrix K th Assembled from the unit heat conduction matrix:

[0089]

[0090] In the formula, ∑ represents the overall matrix assembly operator, is the unit heat conduction matrix. For the finite element problem of a planar four-node unit, the unit heat conduction matrix is ​​defined as:

[0091]

[0092] In the formula, Ω e is the given cell domain, t is the cell thickness, B th is the thermal gradient-temperature matrix, D th is the thermal conductivity matrix.

[0093] Similarly, the global stiffness matrix K is assembled from the element stiffness matrices:

[0094]

[0095] In the formula, ∑ represents the overall matrix assembly operator, k 0 is the unit stiffness matrix. For the finite element problem of a plane four-node unit, the unit stiffness matrix is:

[0096]

[0097] Where B is the strain-displacement matrix.

[0098] Since the thermal expansion caused by the temperature field and the mechanical field load act together on the unit domain, the thermal-mechanical load vector is obtained:

[0099] F=F m +F th

[0100] In the formula, F m is the mechanical load vector, F th is the equivalent thermal load vector. The equivalent thermal load vector can be assembled from the element equivalent thermal load vectors:

[0101]

[0102] Further, from the finite element governing equation:

[0103] KU=F

[0104] The structural displacement field U can be calculated.

[0105] Step 5: Calculate the structural flexibility, global mass fraction, local volume constraint and node temperature constraint values, and perform corresponding sensitivity analysis;

[0106] Specifically, according to the adjoint sensitivity analysis, the flexibility sensitivity expression is:

[0107]

[0108] Where E is the elastic modulus, κ is the thermal conductivity, and f β are thermal stress coefficients, both of which are related to the design variables μ and v.

[0109] In this implementation, the derivative of the global mass constraint with respect to the design variable is:

[0110]

[0111] According to the chain rule, the derivative of the improved mass density ρ′ with respect to the design variable μ is:

[0112]

[0113] Since the improved mass density ρ′ is independent of the design variable v,

[0114]

[0115] The node temperature constraint sensitivity expression is:

[0116]

[0117] In this implementation, since the local volume constraint only acts on the non-uniform filling area, the local volume constraint is only related to the second design variable v, according to the local volume constraint formula:

[0118]

[0119] The local volume constraint sensitivity expression can be obtained as:

[0120]

[0121] Step 6: Determine the objective function and constraints according to the mathematical model of the optimization problem, and submit the sensitivity to the moving asymptotes algorithm (MMA) solver to update the overall topological field μ and the lattice density field v;

[0122] Step 7: Check whether the iteration stop condition is met: the maximum change of the design variable is less than 0.01. If so, output the result; otherwise, return to step 2 and repeat the above steps.

[0123] A specific implementation of the above scheme is as follows:

[0124] This embodiment is to establish a power battery module topology design domain for 18650 cylindrical lithium-ion cells, such as Figure 2 As shown. The rectangular plane size is 150mm×150mm, which is discretized into 300×300 rectangular units. 25 black circular domains are evenly distributed in the plane to represent the location of the 18650 battery cells. The radius of the circular domain is 9mm. All degrees of freedom at the four corners of the design domain are fixed, and uniform loads F are applied to the four boundaries of the design domain. m = 2000N; the concentrated node heat source at the center node of the circular domain is used to replace the battery heat equivalently, and 25 heat source loads P = 1W are applied. The material used is 2018 aluminum alloy, and the material parameters are shown in Figure 3 .

[0125] Step 1: Select the overall topological field μ and the lattice density field v as design variables, establish an optimization model with structural flexibility as the target and mass fraction as the constraint, and initialize it including finite element discretization, definition of mechanical field and temperature field finite element calculation related parameters, PDE filter matrix preprocessing and iteration related parameter initialization;

[0126] Step 2: If Figure 4 As shown, the design variable μ is initialized, and the density field is obtained through two-step PDE filtering and Heaviside projection. and the normalized spatial gradient norm τ, the two filter radii are R1 = 16 and R2 = 10 respectively; the design variable v is initialized, and two-step PDE filtering and Heaviside projection are performed to obtain the non-uniform filling density field ψ, and then the local volume fraction aggregation field is obtained by one-step filtering The two filtering radii are R3=16 and R4=10 respectively; the other parameters such as projection threshold η are as follows Figure 3 shown.

[0127] Step 3: By τ and Substitute the material interpolation function to calculate the elastic modulus E, thermal conductivity κ, and thermal stress coefficient f β and mass density ρ to describe the distribution of material in the design domain;

[0128] Step 4: Assemble the overall thermal conductivity matrix K th And the overall stiffness matrix K, perform thermal-mechanical coupling finite element calculation to solve the structural displacement field U;

[0129] Step 5: Calculate the structural flexibility, global mass fraction, local volume constraint and node temperature constraint values, and perform corresponding sensitivity analysis;

[0130] Step 6: Determine the objective function and constraints according to the mathematical model of the optimization problem, and submit the sensitivity to the moving asymptote algorithm solver to update the overall topological field μ and the lattice density field v;

[0131] Step 7: Check whether the iteration stop condition is met: the maximum change of the design variable is less than 0.01. If so, output the result; otherwise, return to step 2 and repeat the above steps.

[0132] The present invention is not limited to the above-mentioned specific implementation modes. A person skilled in the art may implement the present invention in various other specific implementation modes according to the embodiments and the disclosure of the drawings. Therefore, any design that adopts the design structure and concept of the present invention and makes some simple transformations or changes shall fall within the scope of protection of the present invention.

Claims

1. A thermal-mechanical coupling topology optimization method for a thin-walled non-uniform filling structure, characterized in that: The following steps are involved: Step 1: Establish a finite element model according to the design requirements, determine the non-design domain and design domain, select the overall topology field μ and the lattice density field v as design variables, establish an optimization model with structural flexibility as the target and mass fraction as the constraint, and initialize the model including finite element discretization, definition of mechanical field and temperature field finite element calculation related parameters, PDE filter matrix preprocessing and iteration related parameter initialization; Step 2: Initialize the design variable μ and obtain the density field through a two-step filtering algorithm and the normalized spatial gradient norm τ; initialize the design variable v, perform one-step filtering and projection to obtain the non-uniform filling density field ψ, and then perform one-step filtering to obtain the local volume fraction aggregation field Step 3: By converting the density field Spatial gradient norm τ and local volume fraction aggregation field Substitute the material interpolation function to calculate the elastic modulus E, thermal conductivity κ, and thermal stress coefficient f β and mass density ρ to describe the distribution of material in the design domain; Step 4: Assemble the overall thermal conductivity matrix K th And the overall stiffness matrix K, perform thermal-mechanical coupling finite element calculation to solve the structural displacement field U; Step 5: Calculate the structural flexibility, global mass fraction, local volume constraint and node temperature constraint values, and perform corresponding sensitivity analysis; Step 6: Determine the objective function and constraints according to the mathematical model of the optimization problem, and submit the sensitivity to the moving asymptote algorithm solver to update the overall topological field μ and the lattice density field v; Step 7: Check whether the iteration stop condition is met: the maximum change of the design variable is less than 0.

01. If it is met, output the result; otherwise, return to step 2 and repeat the above steps in sequence; In step 1, the optimization model can be expressed as: s.t.K th T=Q KU=F m +F th lvc(ν)≤0 gmc(μ,ν)=M(μ,v) / M * -1≤0 0≤μ e ≤1,0≤ν e ≤1 In the formula, st represents the constraint condition, U T Where T represents the transpose of the matrix, c is the flexibility, and μ is the relative density of the unit μ e The vector composed of ν is the relative density of the unit ν e The vector composed of K is the overall stiffness matrix of the structure, K th is the overall heat transfer matrix of the structure, U and T are the node displacement vector and node temperature vector respectively, Q is the node heat load vector, F m and F th are the node mechanical load vector and the equivalent node thermal load vector, t i is the temperature value of the ith node, is the temperature constraint value of the ith node, lvc is the local volume constraint of the non-uniform filling area, and gmc is the global mass constraint.

2. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1 is characterized in that: In step 2, the filtering expression is: In the formula, is the Laplace operator, τ is the cell density before filtering, and the cell density after filtering can be obtained by applying the filtering boundary conditions. R f is the filter radius.

3. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1 is characterized in that: In step 2, the cell density after filtering of the e-th cell is The relative density of the unit obtained after Heaviside projection for: Where β is the sharpness of the projection, η is the threshold of the projection, Represents the relative density of the e-th unit after filtering The relative density of the unit obtained after projection; Through the two-step filtering and projection operation on μ, the filling density field is obtained and the thin-wall density field By one-step filtering and projection operation on v, the non-uniform filling density field ψ can be obtained. Applying a smoothing filter on the filling field ψ can obtain the local volume fraction field of the structure:

4. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1, characterized in that: In step 3, the material interpolation function is: Where ρ, E, and k are the overall mass density, elastic modulus, and thermal conductivity, ρ 0 、E 0 and κ 0 They represent the mass density, elastic modulus and thermal conductivity of the shell region material, respectively, ρ , E and λ κ They represent the ratio of the material properties of the filling region to the mass density, elastic modulus and thermal conductivity of the shell region respectively. In order to avoid numerical singularity in the finite element calculation of the empty region material, the minimum quantity ρ is introduced. min 、E min and κ min .

5. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1, characterized in that: The thermal stress coefficient f described in step 3 β The interpolation function is as follows: In the formula, E 0 represents the elastic modulus of the shell region material, λ E It represents the ratio of the material properties of the filling area to the elastic modulus of the shell area material. In order to avoid numerical singularity in the finite element calculation of the empty phase area material, the minimum quantity E is introduced. min .

6. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1, characterized in that: In step 4, the overall thermal conductivity matrix K th And the global stiffness matrix K is assembled as follows: In the formula, ∑ represents the overall matrix assembly operator, k 0 and is the unit stiffness matrix and the unit heat conduction matrix.

7. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1, characterized in that: In step 5, the sensitivity of the flexibility c to the design parameters can be obtained by the following formula: Where E is the elastic modulus, κ is the thermal conductivity, and f β are thermal stress coefficients, both of which are related to the design variables μ and v.

8. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1, characterized in that: In step 5, the sensitivity of the global mass constraint gmc with respect to the design parameters can be obtained by the following formula: According to the chain rule, the derivative of the improved mass density ρ′ with respect to the design variable μ is: Since the improved mass density ρ′ is independent of the design variable v, 9. The thermal-mechanical coupling topology optimization method for thin-walled non-uniform filling structures according to claim 1, characterized in that: In step 5, the node temperature constraint sensitivity and the local volume constraint sensitivity can be obtained by the following formula: x∈{μ,v}

Images