Mechanical-thermal coupling simulation and structural optimization method and system for ultra-complex structures

Abstract

The present invention discloses a method and system for mechanical-thermal coupling simulation and structural optimization of ultra-complex structures. The method comprises: generating parameters and data structures required for subsequent simulation and optimization for a given model and design domain; constructing a background coarse grid, seed points and radii, modeling a porous structure based on a Voronoi structure, and precalculating its stress field; calculating a density field, performing thermodynamic simulation on the generated porous structure based on the density field, the background coarse grid, and a shape function; querying the stress matrix corresponding to each seed point based on the obtained stress field, performing eigendecomposition and extraction on the stress matrix, and obtaining a weight matrix corresponding to the seed point; using a topology optimization method, updating control parameters through a GCMMA method based on sensitivity, thereby updating the porous structure and optimizing material distribution, and continuously optimizing until the change in control parameters is less than a threshold. The method can quickly and effectively solve the topology optimization problem of anisotropic porous structures under multi-physics field conditions.

Description

Technical Field

[0001] This invention belongs to the field of topology optimization technology, and relates to a method and system for force-thermal coupling simulation and structural optimization of ultra-complex structures. Background Technology

[0002] Porous foam materials have attracted widespread attention across various fields due to their unique properties and potential applications. The porous structure of these materials endows them with a wide range of functions, such as enhanced strength-to-weight ratio, thermal insulation, and sound absorption. Therefore, porous foams have become an important component of many industries, including aerospace, automotive, and construction.

[0003] Topology optimization can fully unlock the potential of porous structures, thereby developing structures with superior performance. By systematically redistributing materials within a given design space, topology optimization allows for the identification of optimal configurations that maximize desired properties while satisfying specified constraints. [1,2 Through this process, porous structures can be optimized to achieve improved mechanical, thermal, and acoustic properties.

[0004] However, traditional voxel-based topology optimization methods have both advantages and limitations when applied to porous foam materials. Their advantages include ease of structural modeling and the ability to capture complex geometries. However, these methods are limited in computational efficiency, especially for large-scale optimization problems. [3] Therefore, alternative methods are needed to overcome these limitations, and topology optimization can be used to improve the performance of porous foams.

[0005] Anisotropic structures, where materials exhibit different properties in different directions, hold special significance in various applications. Topology optimization of multiphysics anisotropic porous structures has the potential to promote the development of materials and structures. This optimization can be achieved through various methods, such as thermodynamic orientation optimization. [4] Enhanced hypothetical strain unit [5] and strength-based topology optimization [6] Anisotropy also plays an important role in topology optimization of additive manufacturing load-bearing structures.

[0006] Topology optimization of anisotropic porous structures under multiphysics conditions is currently a hot research topic. (Lind-Nordgren) [7] A method for designing the optimal arrangement of acoustic and vibrational properties of porous open-cell foams using the finite element solution of the Biot equation is discussed. (Harston) [8] A topology optimization method is proposed, which can simultaneously optimize the topology and material properties of anisotropic materials, thereby achieving performance improvement. (Zhao) [9] A novel approach is proposed to optimize the design of porous structures for additive manufacturing and improve their structural performance.

[0007] However, achieving topology optimization of anisotropic porous structures under multiphysics conditions remains a significant challenge. Its complexity lies in considering both various physical factors and the anisotropic design of the porous structure. Achieving optimization that takes into account the anisotropic behavior of materials in different directions adds further complexity. The interactions between different physical fields further complicate the optimization process, necessitating the development of innovative methods to effectively address these challenges.

[0008] [1]J.Zhu,J.Sun,H.Tang,J.Wang,Q.Ao,T.Bao,W.Song,Gradient-structuraloptimization of metal fiber porous materials for sound absorption,PowderTechnology 301(2016)1235–1241.doi:10.1016 / j.powtec.2016.08.006.

[0009] [2] J. Yoon, H. Kim, T. Koh, S. Pyo, Microstructural characteristics of soundabsorbable porous cement-based materials by incorporating natural fibers and aluminum powder, Construction and Building Materials 243(2020)118167.doi:10.1016 / j.conbuildmat.2020.118167.

[0010] [3] S.Dou, Aprojection approach for topology optimization of porousstructures through implicit local volume control, Structural and Multidisciplinary Optimization 62(2)(2020)835–850.doi:10.1007 / s00158-020-02539-x.

[0011] [4]D.R.Jantos,K.Hackl,P.Junker,Topology optimization with anisotropicmaterials,including a filter to smooth fiber pathways,Structural andMultidisciplinary Optimization 61(5)(2020)2135–2154.doi:10.1007 / s00158-019-02461-x.

[0012] [5]G.Zhang,L.Li,K.Khandelwal,Topology optimization of structures withanisotropic plastic materials using enhanced assumed strain elements,Structural and Multidisciplinary Optimization 55(6)(2016)1965–1988.doi:10.1007 / s00158-016-1612-1.

[0013] [6]A.M.Mirzendehdel,B.Rankouhi,K.Suresh,Strength-based topologyoptimization for anisotropic parts,Additive Manufacturing 19(2018)104–113.doi:10.1016 / j.addma.2017.11.007.

[0014] [7]E.Lind-Nordgren,P. Optimising open porous foam foracoustical and vibrational performance,Journal of Sound and Vibration 329(7)(2010)753–767.doi:10.1016 / j.jsv.2009.10.009.

[0015] [8] S.Harston, C.Mattson, M.Koecher, Atopology optimization method with anisotropic materials, in: 13th AIAA / ISSMO Multidisciplinary Analysis Optimization Conference, American Institute of Aeronautics and Astronautics, 2010.doi:10.2514 / 6.2010-9176.

[0016] [9] Summary of the Invention

[0017] To address the shortcomings of existing technologies, this invention provides a method and system for mechanical-thermal coupling simulation and structural optimization of ultra-complex structures, aiming to solve the topology optimization problem of anisotropic porous structures under multiphysics conditions. Specifically, the contributions of this invention can be summarized as follows:

[0018] 1) A porous structure modeling method with control parameters of clear design significance is provided to accurately and robustly simulate porous foam.

[0019] 2) A simulation method based on a fixed background mesh and CBN shape function is provided to achieve efficient thermodynamic simulation of the model.

[0020] 3) A comprehensive topology optimization framework based on control parameters is provided, which is suitable for the design optimization of multi-objective functions of porous foam.

[0021] 4) An anisotropic porous foam optimization method was implemented. By guiding the anisotropic optimization of the porous structure with a fixed density field, the mechanical properties and regularity of the final porous foam were optimized.

[0022] Based on the above contributions 1)-4), this invention presents a design scheme for topology optimization of anisotropic porous structures under multiphysics conditions.

[0023] The technical solution adopted in the present invention is as follows:

[0024] A method for force-thermal coupling simulation and structural optimization of ultra-complex structures includes:

[0025] 1) For the given model and design domain, generate and establish the parameters and data structures required for subsequent simulation and optimization; including constructing a background coarse mesh, constructing seed points and radii, modeling porous structures based on the Voronoi structure, and pre-calculating their stress fields.

[0026] 2) Calculate the density field. Based on the obtained density field, perform thermodynamic simulation on the generated porous structure using the background coarse mesh and shape function. During the thermodynamic simulation, calculate the parameter values ​​and sensitivity of the relevant objective function and constraint conditions.

[0027] Based on the obtained stress field, the stress matrix corresponding to each seed point is queried, and the stress matrix is ​​subjected to feature decomposition and extraction processing to obtain the weight matrix corresponding to the seed point. The objective function value and sensitivity characterizing the anisotropy of the porous structure are calculated.

[0028] Using topology optimization, control parameters are updated based on sensitivity using the GCMMA method, thereby updating the porous structure and optimizing material distribution;

[0029] 3) Repeat step 2) until the change in control parameters is less than the threshold. After the above process is completed, the final seed point and radius, as well as the topological information of the porous structure, are obtained. The final optimized porous foam is constructed based on the directed distance field using the above structure.

[0030] A force-thermal coupling simulation and structural optimization system for ultra-complex structures is provided to implement the method described above.

[0031] The beneficial effects of this invention are:

[0032] 1) The parameterized model generated by the porous structure modeling method provided by this invention can accurately simulate porous foam. Compared with the voxel model, it improves the simulation and optimization efficiency while ensuring the characteristic parameters.

[0033] 2) The simulation method based on a fixed background mesh and CBN shape function used in this invention is more efficient than traditional mesh-based simulation methods.

[0034] 3) The topology optimization framework proposed in this invention avoids expensive global sensitivity calculations and significantly improves computational efficiency.

[0035] 4) The anisotropic porous foam optimization method implemented in this invention has better mechanical properties compared with the isotropic structure optimization method. Attached Figure Description

[0036] Figure 1 This is a schematic diagram of the overall framework of the method of the present invention;

[0037] Figure 2 The definition of beams in implicit modeling of porous structures includes descriptions of multiple beams and beam radii, and descriptions of single beams.

[0038] Figure 3 This is a schematic diagram of the design domain and boundary in a specific embodiment of the present invention;

[0039] Figure 4 This is a porous foam generated based on the final optimization result in a specific embodiment of the present invention;

[0040] Figure 5 This is a simulation of the stress field distribution of the finally optimized porous structure in a specific example of the present invention;

[0041] Figure 6 This is a temperature field distribution obtained from the simulation of the finally optimized porous structure in a specific embodiment of the present invention;

[0042] Figure 7 This is a convergence curve of the objective function and constraints in a specific embodiment of the present invention;

[0043] Figure 8 A comparison of flexibility (stiffness) under anisotropic and isotropic conditions in a specific embodiment of the present invention;

[0044] Figure 9 This is a comparison of the average temperature (heat dissipation) under anisotropic and isotropic conditions in a specific embodiment of the present invention. Detailed Implementation Plan

[0045] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0046] According to a specific embodiment of the present invention, such as Figure 1 As shown, a method for force-thermal coupling simulation and structural optimization of ultra-complex structures includes the following:

[0047] First, given the model and design domain, the parameters and data structures required for subsequent simulation and optimization are generated and established. For example, a coarse background mesh is constructed through tetrahedral meshing, and seed points for constructing porous structures are generated and initialized. Simultaneously, a voxel-based topology optimization method is used to pre-calculate the stress field of the model, obtaining the stress field that guides anisotropic optimization.

[0048] By constructing a background coarse mesh, seed points, and radii, a Voronoi structure is generated within the design domain, thus modeling a porous structure composed of multiple beams with thicknesses. Further, the distance from each mesh element node to the nearest beam is calculated, and based on the density function, this distance is transformed into the density of each node. The average density of all nodes in that mesh element is then calculated to obtain the density field of all mesh elements.

[0049] Based on the obtained density field, thermodynamic simulations were sequentially performed on the generated porous structure using a simulation method based on the background mesh and CBN shape function. Based on the temperature field obtained from the thermal simulation, the force loads exerted on the porous structure under thermal conditions were calculated. Independent force loads applied under combined force conditions were then applied to the porous structure, followed by mechanical simulations. During the thermodynamic simulation, the parameter values ​​and sensitivities of the relevant objective function and constraint conditions were calculated.

[0050] Based on the obtained stress field, the stress matrix corresponding to each seed point is queried, and the stress matrix is ​​subjected to feature decomposition and extraction to obtain the weight matrix corresponding to the seed point. The objective function value and sensitivity characterizing the anisotropy of the porous structure are then calculated.

[0051] Furthermore, using topology optimization methods, the control parameters are updated based on the sensitivity using the GCMMA method, thereby updating the porous structure and optimizing the material distribution.

[0052] Repeat the above steps until the control parameters change very little or not at all.

[0053] Once the above process is completed, the final seed point and radius, as well as the topological information of the porous structure, will be obtained. Using this structure, the final optimized porous foam is constructed based on a directed distance field.

[0054] In the above scheme, in the modeling part, this invention abstracts the model into a beam structure network controlled by the position and radius of seed points. To achieve explicit shape control, this invention utilizes the edges of the Voronoi diagram with specific radii to form Voronoi foam.

[0055] In the simulation section, to avoid complex mesh generation, this invention introduces a fixed background mesh. This invention explores efficient multiphysics simulation under force-thermal coupling conditions and temperature constraints.

[0056] In the optimization section, this invention studies the generation of anisotropic porous structures in order to maintain a certain regularity while preserving physical properties, thereby improving the manufacturability of the model.

[0057] Based on this, this invention establishes a topology optimization framework that can achieve anisotropic mechanical-thermal coupling topology optimization of porous structures based on Voronoi generation and implicit representation. It is solved according to the classical topology optimization framework, and follows the... Figure 1 The method described is as follows. This method first performs downflow simulation using a coarse background mesh embedded in the domain, and then pre-calculates the stress field using a voxel topology optimization method. After initializing the design variables X and radius r, the variable values ​​are repeatedly updated according to subsequent steps until convergence: density function calculation, property simulation, sensitivity calculation, and design update. Finally, the optimized Voronoi porous structure is constructed using X, r, Ω, and the output.

[0058] In density function calculations, this invention constructs an implicit function to represent the Voronoi foam structure and generates smooth foam under explicit geometric control variables. The implicit representation avoids the unreliable and time-consuming geometric calculations involved in explicitly constructing Voronoi foam and generating simulation meshes, and allows for high-precision continuous finite element integration.

[0059] In property simulation, this invention introduces material-aware CBN shape functions to reduce computational costs. Simulations were performed on the initially constructed embedded coarse mesh, and high simulation accuracy was still achieved by trimming a set of polyhedral material-aware shape functions for each coarse element. Furthermore, this process directly operates on the implicit foam representation, thereby eliminating the need for error-prone boundary conformal fine mesh generation. This method significantly reduces computational costs.

[0060] Furthermore, to achieve anisotropic Voronoi foam optimization, this invention utilizes the weight matrix Mx in the pre-calculated objective function of the stress field to achieve field-guided anisotropic optimization. Compared to traditional isotropic canonical optimization (when Mx is an identity matrix), the method of this invention can obtain anisotropic Voronoi cells, thereby improving mechanical properties while maintaining regularity.

[0061] The technical details used in the solution are explained in detail below.

[0062] 1. Implicit Modeling of Porous Foams

[0063] First, an open Voronoi foam is defined as the edge (with a certain radius) of a Voronoi cell confined within the model. By mixing together the implicit functions of all beams introduced by the Voronoi edges, or together with the shell, an implicit representation of the Voronoi foam defined as a density field can be obtained. Details are as follows.

[0064] Each Voronoi foam Ω V(X,r) consists of two parts: an internal part and a boundary part. For the internal part, this invention first starts from the seed X = {X...} i i = 1, ..., N s Calculate the Voronoi diagram. Within a sufficiently large bounding box, cut all its edges using the outer free shape Ω. The remaining j-th edge becomes a beam with a radius equal to the average of its neighboring seed points.

[0065]

[0066] Where χ j Let |χ be the set of adjacent seeds of edge j. j | represents the number of seeds. The radius of the seeds is collected in the form of a vector r, such as... Figure 2 As shown in the left figure.

[0067] A beam φ is defined corresponding to the Voronoi edges of vertices v1 and v2, consisting of a radius of... A cylinder of height ||v1-v2|| and two radii of Composition of the hemispherical end: see Figure 2 As shown on the right. Implicit form of beam φ. for:

[0068]

[0069] In the formula, d(x,v2,v2) represents the minimum distance from point x to edges v1 and v2.

[0070]

[0071] in

[0072]

[0073] The implicit representation of the entire open-cell Voronoi foam is achieved by merging the implicit functions of all beams after smoothing using the Kreisselmeier-Steinhauser (KS) function.

[0074]

[0075] In the formula φ max =max(φ1,…,φ n ), where n is the number of beams, and p = 16. The KS function makes the describing function Φ(x) compact and differentiable with respect to x.

[0076] The boundary portion associated with the external shape Ω is also implicitly represented by φ. Ω (x) represents this. When a closed shell is required, φ is... Ω(x) can be directly extended to equation (5) as an additional beam. Otherwise, the invention includes the intersection between the outer shape and the Voronoi surface, where the thickness of each surface is the average of the two seed radii that determine the surface. For simplicity, the invention will not distinguish between shelled and shellless Voronoi foams.

[0077] Voronoi Foam Ω V (X,r) is ultimately described as a density field by the regularized Heaviside function H(Φ(x)).

[0078]

[0079] Where ∈ controls the magnitude of regularization, and α = 1e -6 The default value is a small positive number to avoid singularity in the global stiffness matrix.

[0080] 2. Finite element formula

[0081] To achieve the simulation of implicit foam, this invention immerses the model in a fixed tetrahedral background mesh D = {D e The shape Ω is covered in the sequence {e = 1, 2, ..., N}, i.e. This not only avoids the time-consuming and error-prone model partitioning process, but also helps to achieve convergence in subsequent optimizations by fixing the mesh.

[0082] In each tetrahedral element D e In the middle, density H e Let H(x) be the average value of its four nodes.

[0083]

[0084] Where, For unit D e The coordinates of the k-th node.

[0085] This invention considers a linear elastic problem involving coupled heat conduction and thermal expansion effects. It considers a linear elastic problem with a fixed displacement boundary Γ. u and fixed temperature boundary Γ T The design domain Ω. For example... Figure 3 As shown, surface forces act on the Neumann boundary Γ N Above, the surface heat flux is located at the boundary Γ h Above. The fundamental equation for the force-thermal coupling problem is:

[0086]

[0087] In the formula, u is the displacement vector, σ is the stress vector, T is the temperature, and T0 is the Γ. T The temperature specified above, ΔT = T – T ref T refWhere is the thermal expansion reference temperature, k is the thermal conductivity, f is the load vector, and Q is the internal strength of the heat source.

[0088] The stress vector is represented as:

[0089] σ=D(ε-ε th (9)

[0090] In the formula, D is the elastic matrix, ε is the strain vector, and ε th This is the thermal strain vector.

[0091] The strain energy E can be expressed as:

[0092]

[0093] After discretizing the design domain using the finite element method, the displacement vector u and the stress tensor σ can be expressed as the relationship between the nodal displacements U:

[0094]

[0095] In the formula, B is the strain-displacement matrix, and N is the shape function matrix. E is rewritten as:

[0096]

[0097] The total potential energy can be expressed as

[0098]

[0099] In the formula f v f is the force vector of the object. b For boundary traction force.

[0100] According to the principle of minimum potential energy, the change in the nodal displacement vector U can be expressed as:

[0101]

[0102] Consider the variational principle, when f v When =0, the global system of equations can be obtained from Equation 15:

[0103] K m U = F = F m + F th (16)

[0104] have

[0105]

[0106] Where K m Let F be the stiffness matrix, U be the nodal displacement vector, and F be the stiffness matrix. m F is the mechanical load vector. th This is the thermal load vector.

[0107] The static equilibrium equations of a linear thermoelastic continuum under mechanical and temperature loads can be expressed as Equations 16 and 17. m From the element stiffness matrix K m,e Combination

[0108] K m =∑ e K m,e (18)

[0109] K m,e It can be calculated as an integral.

[0110]

[0111] In the formula B e Ω is the element strain-displacement matrix. e This is a unit physical domain. Similarly, the thermal load F th The unit thermal load F th,e Combination

[0112] F th =∑ e F th,e (20)

[0113] in

[0114]

[0115] In the formula ε th,e Let T be the element's thermal strain vector, α be the coefficient of thermal expansion, and T be the thermal strain vector. e The element temperature is T, which is the element node temperature. e The average value of Φ is a vector [1 1 1 0 0 0] in three-dimensional space. Similarly, the displacement u(x) in the background mesh D can be obtained by interpolating the nodal displacements U.

[0116] u(x) = N e (x)U e (twenty two)

[0117] In the formula N e (x) is the shape function matrix at the unit node, U e This represents the displacement of the element nodes.

[0118] In the case of force-thermal coupling, the effect of force on heat is minimal. To simplify the analysis, this invention focuses on the unidirectional coupling of heat and force. Specifically, this invention incorporates the thermal stress generated in the material into the mechanical equations.

[0119] The heat balance equation for steady-state heat conduction is:

[0120] KT (X,r)T=P (23) where K T Let K be the thermal conductivity matrix, T be the nodal temperature vector, and P be the heat flux density vector. T And P is determined by the unit thermal conductivity matrix K T,e and unit heat flux vector P e Combination

[0121]

[0122] The thermal conductivity matrix can be derived using the Galerkin finite element method.

[0123]

[0124] In the formula λ e The effective thermal conductivity.

[0125] By interpolating the nodal temperature T, the temperature T(x) in the background mesh D can be obtained.

[0126] T(x) = N e (x)T e (26)

[0127] Where T e It is the unit node temperature, which satisfies

[0128]

[0129] 3. Optimize problem design

[0130] The Voronoi foam design problem uses the seed and radius X,r as design variables for performance optimization. Its analytical and discrete interpretations are as follows.

[0131] This invention aims to create Voronoi foam Ω V The overall deformation of (x,r) is small, while maintaining a certain shape regularity. Therefore, this invention introduces the terms "compliance" and "shape". The compliance C(X,r,U), widely used in topology optimization to measure the elastic potential of an object, is set as the physical objective.

[0132] C(X,r,U)=F T U (28)

[0133] The shape regulation energy S(X) is used to regularize Voronoi cells into approximately regular polyhedra.

[0134]

[0135] That is, each Voronoi cell Seed X i With center of mass The Euclidean distance between them, via The weighted sum.

[0136] To obtain the anisotropic structure of Voronoi foam, this invention introduces M... x As a weight matrix for any coordinate x, different weights are assigned to different directions and Voronoi cells. The following section will introduce M in detail. x The calculation process is as follows: First, simulation and optimization are performed based on a high-resolution voxel model to obtain the stress field of the optimized model. The Cauchy stress tensor at point x in the stress field is obtained from...

[0137]

[0138] Where σ ii For normal stress, σ ij ,i≠j represents shear stress. Due to the symmetry of shear stress, σ ij =σ ji Therefore σ x It is a symmetric matrix. To extract σ x The principal direction vector is decomposed into its eigenvalues ​​in this invention:

[0139] σ x =VΛV -1 (31)

[0140] Where Λ is a diagonal matrix with its eigenvalues ​​on the diagonal, and V is a matrix with its eigenvectors as its columns. This invention extracts the eigenvector corresponding to the largest eigenvalue as the principal direction vector v. m .

[0141] v m =[v1 v2 v3] (32)

[0142] M x By v m Components:

[0143]

[0144] For simplicity, for each Voronoi unit Take M x x is the center of mass. As shown in Equation 29.

[0145] Finally, for the optimization problem with volume and temperature constraints, the optimal design parameters p = {p} are obtained. j Let p = {j = 1, 2, ..., n}, where p = ... j Represents a component of variable X or r:

[0146]

[0147] In this example, the objective function J(p,U) is J(X,r,U), which is set as the weighted sum of the physical objective and the shape regularization term.

[0148] J(p, U) = (1-w) C(p, u) + w S(X), (35)

[0149] V and For the total volume and maximum volume constraints of the structure, T i and For the node temperature and the corresponding maximum temperature constraint, p j and Design variable p j The lower and upper bounds.

[0150] As mentioned above, different physical properties or shape regularization measures can be introduced for different design purposes, and the target foam can be derived similarly in this process.

[0151] To solve the optimization problem Eq.(34), the Global Converging Moving Asymptote Method (GCMMA) based on numerical gradient is adopted due to its robust convergence in design optimization. The original non-convex problem is approximated by a set of convex subproblems by utilizing the gradients of the optimization objective and constraints with respect to the design variable p.

[0152] 4. Sensitivity Analysis

[0153] 1.1 Sensitivity to Flexibility

[0154] To optimize the objective function, we first consider the effect of compliance on the design parameter p. j Sensitivity:

[0155]

[0156] In order to eliminate Eq.(36) The partial derivative of this invention with respect to Eq.(23) is:

[0157]

[0158] Transposing these terms gives us:

[0159]

[0160] Similarly, taking the partial derivative with respect to Eq.(16) and rearranging the terms:

[0161]

[0162] For simplicity, assume P and p j Irrelevant, let's assume Substitute Eq.(38) and Eq.(39) into Eq.(36).

[0163]

[0164] To simplify Eq.(40), λ t Defined as neighbor vectors, which satisfy:

[0165]

[0166] λ t This can be obtained by solving the following system of linear equations:

[0167]

[0168] According to Eq.(16) and Eq.(21),

[0169]

[0170] λ t Substituting into Eq.(40), we finally obtain the compliance C relative to the design parameter p. j Sensitivity:

[0171]

[0172] According to Eq.(25), Eq.(19) and Eq.(21), the present invention has

[0173]

[0174]

[0175] and

[0176]

[0177] in

[0178]

[0179] 1.2 Temperature Constraint Sensitivity

[0180] To meet temperature constraints, this invention will next calculate the sensitivity of temperature relative to design parameters. Node temperatures are written as...

[0181]

[0182] in This represents a vector where all values ​​except the i-th element are 0. Node temperature affects design parameter p. j The gradient is expressed as

[0183]

[0184] Differentiate Eq.(23)

[0185]

[0186] Substituting Eq.(51) into Eq.(50), we get

[0187]

[0188] To simplify Eq.(52), similar to λ in Eq.(41) t This invention introduces an adjacent vector λ i This can be obtained by solving the following system of linear equations:

[0189] K T λ i =-L i (53)

[0190] Then, the sensitivity to node temperature is calculated.

[0191]

[0192] 5. Efficient calculation of CBN shape function

[0193] Traditional finite element methods require very fine mesh elements to achieve high-precision simulations, which introduces a large number of degrees of freedom. This poses a challenge to solving linear equations.

[0194] Building upon previous research, this invention introduces a material-aware CBN shape function to replace the original shape function, thereby reducing the degrees of freedom and accelerating computation.

[0195] First, the entire virtual domain is coarsely meshed. For each coarsened unit e and its virtual domain Ω... e The present invention further divides it into a tetrahedral mesh, obtaining a refined mesh h of e. To reduce the computational degrees of freedom while maintaining computational accuracy, the present invention uses [a specific method] to calculate simulation values ​​such as the stiffness matrix. Subsequently, the present invention compresses the computational results by mapping the calculation results to a coarse mesh.

[0196] The unit shape function N in Eq.(26) and Eq.(22) e By transforming the matrix Replace with the following cell shape functions, which consist of linear shape functions on a fine mesh:

[0197]

[0198] Where N e,h (x) is Ω eThe set of node shape functions on a fine mesh. Similarly, the strain-displacement matrix B in Eq.(25), Eq.(19), and Eq.(21) e It can be replaced with

[0199]

[0200] In the formula B e,h Ω e Global strain-displacement matrix for fine mesh. CBN transformation matrix. The aim is to map CBN node values ​​to Ω. e The internal values ​​are derived from the boundary interpolation matrix Ψ and the boundary-interior transformation matrix. The product is obtained as follows:

[0201]

[0202] Among them Ψ and Map the displacements of CBN to the boundary nodes, and then to Ω. e Fully fine-grained nodes in the array.

[0203] The boundary interpolation matrix Ψ maps CBN node values ​​to fine-grid boundary node values ​​Ω. e It is derived by constructing a bicubic B′ezier interpolation surface on the surface of interest, using CBN as control points. The matrix Ψ is formed by calculating the surfaces at the fine mesh nodes within the surface and collecting them row by row into a matrix, representing all values ​​of the mesh h boundary.

[0204] The transformation matrix maps the boundary node values ​​to Ω. e The node values ​​of a medium-fine mesh. It is calculated using the equilibrium equations for Ω. e Obtained through local simulation using fine mesh

[0205]

[0206] Where k b ,k i ,k bi ,k ib Ω e Upper local stiffness matrix k e submatrix, q b ,q i The vectors for the boundary and internal nodes are respectively, f b This represents the exposed force vector at the boundary node formed by harmonic analysis.

[0207] From the second line, the present invention obtains q i =M e q b The relation, in which Therefore, qi and q b The set is q = [q b ,q i ] T The present invention has form

[0208]

[0209] Among them I 2b It is a 2b×2b identity matrix.

[0210] Therefore, the original coarse element thermal conductivity matrix K T,e Given by Eq.(25), using Replace B e ,Right now

[0211]

[0212] Among them, K T,e,h For fine grid Ω e High-fidelity thermal conductivity matrix, H e,α Ω is the average value of H(x) on tetrahedral elements in the fine mesh h. e,α Let be the imaginary field of sub-element α in element e. Similarly, in the coarsening simulation, the gradient of the element thermal conductivity matrix in Eq.(45) is calculated as follows:

[0213]

[0214] Similarly, K m,e ,F th,e Its gradient can be replaced with:

[0215]

[0216] and

[0217]

[0218] In the formula K m,e,h and F th,e,h Ω e High-fidelity stiffness matrix and thermal load of fine mesh.

[0219] To verify the effectiveness of the method proposed in this invention, experimental results are presented and explained below.

[0220] Under conditions of volume and temperature constraints, topology optimization experiments with minimum flexibility are conducted. Figure 4 , 5Figures 6 and 7 respectively demonstrate the optimal design, stress, and temperature distribution of the optimized porous structure. It can be seen that a clear topological structure was obtained through optimization, and the temperature requirements of the heat source were met. Furthermore, thanks to the temperature constraint of the heat source, the maximum temperature of the structure can be controlled within the permissible value of 325℃. Figure 7 The iterative process of the topology optimization problem is illustrated, demonstrating that the method of this invention can achieve stable convergence. Specifically, the change in the objective value compliance and the change in the material volume fraction exhibit a synchronous fluctuation pattern (see reference). Figure 7 (Comparison of iteration history) aligns with the traditional minimum compliance optimization case based on volume constraints. Furthermore, the method provided by this invention is highly efficient, requiring only 200 iterations to obtain the optimal design. According to the given iteration history, convergence was actually achieved on the 50th iteration.

[0221] Subsequently, the obtained optimal topology was compared with the optimal topology obtained under isotropic conditions. The relevant data for the optimization results under isotropic and anisotropic conditions are as follows: Figure 8 and 9 As shown. Figure 8 The iterative history of compliance under two conditions is presented, demonstrating that the structure generated considering anisotropy possesses superior physical properties. The iterative changes in average temperature during the optimization process (e.g., Figure 9 The study found that the structure generated considering anisotropy has a lower average temperature, which means that the structure optimized considering anisotropy has better heat dissipation performance.

[0222] In summary, the method proposed in this invention demonstrates significant advantages in improving design efficiency, ensuring constraint satisfaction, optimizing structural performance, and adapting to complex material properties, providing a powerful tool for design optimization in related fields.

Claims

1. A method for force-thermal coupling simulation and structural optimization of ultra-complex structures, characterized in that, include: 1) Given the model and design domain, generate and establish the parameters and data structures required for subsequent simulation and optimization; This includes constructing a background coarse mesh, constructing seed points and radii, modeling porous structures based on the Voronoi structure, and pre-calculating their stress fields; 2) Calculate the density field, and perform thermodynamic simulation on the generated porous structure based on the obtained density field, background coarse mesh, and shape function; During the thermodynamic simulation, the parameter values ​​and sensitivities of the relevant objective function and constraint conditions are calculated. Based on the obtained stress field, the stress matrix corresponding to each seed point is queried, and the stress matrix is ​​subjected to feature decomposition and extraction processing to obtain the weight matrix corresponding to the seed point. The objective function value and sensitivity characterizing the anisotropy of the porous structure are calculated. Using topology optimization, control parameters are updated based on sensitivity using the GCMMA method, thereby updating the porous structure and optimizing material distribution; 3) Repeat step 2) until the change in control parameters is less than the threshold. After the above process is completed, the final seed point and radius, as well as the topological information of the porous structure, are obtained. The final optimized porous foam is constructed based on the directed distance field using the above structure.

2. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 1, characterized in that, In step 1), a background coarse mesh is constructed by meshing based on tetrahedrons. The stress field is pre-calculated using a topology optimization method based on voxel method to obtain the stress field that guides anisotropic optimization.

3. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 1, characterized in that, In step 1), the model is abstracted into a beam structure network controlled by the position and radius of the seed points. Seed points for constructing the porous structure are generated and initialized. By constructing the background coarse mesh, seed points and radii, Voronoi foam structure is generated in the design domain, thereby modeling a porous structure composed of multiple beams with thickness.

4. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 3, characterized in that, Implicit functions are constructed to represent Voronoi foam structures. Open Voronoi foam is defined as Voronoi cells with edges of a certain radius confined within the model. The implicit functions of all beams introduced by the Voronoi edges are mixed together, or combined with the shell, to obtain an implicit representation of Voronoi foam defined as a density field. Specifically: Each Voronoi foam Ω V (X,r) consists of two parts: the interior part and the boundary part; For the internal portion, within the bounding box, all its edges are cut using the external freeform shape. The remaining j-th edge becomes a beam with radius equal to the average of the adjacent seed points of edge j, corresponding to the Voronoi edges of vertices v1 and v2. The beam φ is defined by a radius of... A cylinder of height ||v1-v2|| and two radii of The hemispherical ends are composed of the implicit form of beam φ. for: In the formula, d(x,v1,v2) represents the minimum distance from point x to edges v1 and v2. in a=v2-v1,b=x-v1,e=x-v2, The implicit representation of the entire open-cell Voronoi foam is achieved by merging the implicit functions of all beams after smoothing using the Kreisselmeier-Steinhauser (KS) function. In the formula φ max =max(φ1,…,φ n ), where n is the number of beams, and p is 16; The boundary portion associated with the external shape Ω is also implicitly represented by φ. Ω (x) represents; Voronoi Foam Ω V (X,r) is ultimately described by the regularized Heaviside function H(Φ(x)) as a density field: Where ∈ controls the magnitude of regularization, and α = 1e -6 .

5. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 4, characterized in that, In step 2), to simulate implicit foam, the model is immersed in a fixed tetrahedral background mesh D = {D}. e The shape Ω is covered in the sequence {e = 1, 2, ..., N}, i.e. In each tetrahedral element D e In the middle, density H e Let H(x) be the average value of its four nodes. In the formula, For unit D e The coordinates of the k-th node; Consider a boundary Γ with a fixed displacement. u and fixed temperature boundary Γ T The design domain Ω, surface forces act on the Neumann boundary Γ N Above, the surface heat flux is located at the boundary Γ h superior; The static equilibrium equation of a linear thermoelastic continuum under mechanical and temperature loads is expressed as: K m U=F=F m +F th have K m =∫ Ω B T DBdΩ F th H∫ Ω B T Dε th dΩ Where K m Let F be the stiffness matrix, U be the nodal displacement vector, and F be the stiffness matrix. m F is the mechanical load vector. th This is the thermal load vector; B is the strain-displacement matrix, N is the shape function matrix; ε th This is the thermal strain vector; K m From the element stiffness matrix K m,e Combination K m =∑ e K m,e In the formula B e Ω is the element strain-displacement matrix. e For a single physical domain; thermal load F th The unit thermal load F th,e Combination F th =∑ e F th,e e th =α(T e -T ref )F T In the formula ε th,e Let T be the element's thermal strain vector, α be the coefficient of thermal expansion, and T be the thermal strain vector. e The element temperature is T, which is the element node temperature. e The average value, T ref Let Φ be the thermal expansion reference temperature, and Φ be a vector in three-dimensional space [1 1 1 0 0 0]. Similarly, the displacement u(x) in the background mesh D is obtained by interpolating the nodal displacement U. u(x)=N e (coin e In the formula N e (x) is the shape function matrix at the unit node, U e For element node displacement; Incorporating the thermal stress generated in the material into the mechanical equations, the heat balance equation for steady-state heat conduction is K. T (X,r)T=P In the formula K T Let T be the thermal conductivity matrix, P be the nodal temperature vector, and K be the heat flux density vector. T And P is determined by the unit thermal conductivity matrix K T,e and unit heat flux vector P e Composed of various elements, including K T (X,r)=∑ e K T,e (X,r) P=∑ e P e In the formula λ e For effective thermal conductivity, The temperature T(x) in the background mesh D is obtained by interpolating the nodal temperature T. T(x)=N e (x)T e 。 6. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 5, characterized in that, The Voronoi foam design problem uses seed X and radius r as design variables for performance optimization, with the physical objective being the compliance C(X,r,U) of the elastic potential of the object. C(X,r,U)=F T U Voronoi cells are regularized into approximately regular polyhedra by using shape-modulating energy S(X). That is, each Voronoi cell Seed X i With center of mass The Euclidean distance between them, expressed by the weight matrix Weighted sum; For an optimization problem with volume and temperature constraints, find the optimal design parameters p = {p j Let p = {j = 1, 2, ..., n}, where p = ... j Represents a component of variable X or r: The objective function J(p,U) is set as the weighted sum of the physical objective and the shape regularization term. J(p,U)=(1-w)C(p,u)+w S(X), V and For the total volume and maximum volume constraints of the structure, T i and For the node temperature and the corresponding maximum temperature constraint, p j and Design variable p j The lower and upper bounds; The optimization problem is solved by using the gradients of the optimization objective and constraints with respect to the design variable p, employing the Global Converging Moving Asymptote Method (GCMMA) based on numerical gradients.

7. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 6, characterized in that, Flexibility C versus design parameter p j The gradient is: have in The temperature gradient relative to the design parameters is:

8. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 6, characterized in that, For each Voronoi unit Where the weight matrix at point x The calculation is as follows: Based on the stress field obtained through voxel model simulation and optimization, the Cauchy stress tensor at point x in the stress field has: Where σ ii For normal stress, σ ij ,i≠j represents shear stress. Due to the symmetry of shear stress, σ ij =σ ji Therefore σ x Given a symmetric matrix, decompose its eigenvalues: s x =VΛV -1 Where Λ is a diagonal matrix with its eigenvalues ​​on the diagonal, and V is a matrix whose eigenvectors are its columns. The eigenvector corresponding to the largest eigenvalue is extracted as the principal direction vector v. m : v m =[v1v2v3] M x By v m Components:

9. The method for force-thermal coupling simulation and structural optimization of ultra-complex structures according to claim 6, characterized in that, Material-aware CBN shape functions are introduced to replace the original shape functions, thereby reducing degrees of freedom and accelerating computation. Specifically, for the entire virtual domain divided by a coarse mesh, each coarsened element e and its virtual domain Ω are... e Furthermore, it is divided into a tetrahedral mesh to obtain the refined mesh h of e; The unit shape function N e Transformation matrix via CBN Replace with cell shape functions consisting of linear shape functions on a fine mesh: Where N e,h (x) is Ω e A set of node shape functions on a fine mesh; Strain displacement matrix B e Replace with In the formula B e,h Ω e Global strain-displacement matrix of fine mesh; CBN transformation matrix The aim is to map CBN node values ​​to Ω. e The internal values ​​are determined by the boundary interpolation matrix Ψ and the boundary-interior transformation matrix. The product is obtained as follows: Among them Ψ and Map the displacements of CBN to the boundary nodes, and then to Ω. e The fully fine-grained nodes in the model; similarly, in the simulation, the gradient of the element thermal conductivity matrix is ​​calculated as follows: Similarly, K m,e ,F th,e Its gradient is replaced with: and In the formula K m,e,h and F th,e,h Ω e High-fidelity stiffness matrix and thermal load of fine mesh.

10. A system for force-thermal coupling simulation and structural optimization of ultra-complex structures, characterized in that, Used to implement the method as described in any one of claims 1-9.

Images