CN117373579B - Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain - Google Patents

Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain Download PDF

Info

Publication number
CN117373579B
CN117373579B CN202311492300.1A CN202311492300A CN117373579B CN 117373579 B CN117373579 B CN 117373579B CN 202311492300 A CN202311492300 A CN 202311492300A CN 117373579 B CN117373579 B CN 117373579B
Authority
CN
China
Prior art keywords
microstructure
scale
macroscopic
time domain
design
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202311492300.1A
Other languages
Chinese (zh)
Other versions
CN117373579A (en
Inventor
江旭东
张巍
武子旺
滕晓艳
马佳琪
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin University of Science and Technology
Original Assignee
Harbin University of Science and Technology
Filing date
Publication date
Application filed by Harbin University of Science and Technology filed Critical Harbin University of Science and Technology
Priority to CN202311492300.1A priority Critical patent/CN117373579B/en
Publication of CN117373579A publication Critical patent/CN117373579A/en
Application granted granted Critical
Publication of CN117373579B publication Critical patent/CN117373579B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Abstract

A multi-microstructure multi-scale parallel topology optimization method and system under the action of time domain dynamic load belong to the field of structure optimization design, a double Helmholtz smooth-blocking projection scheme is introduced in the method, macroscopic equivalent mechanical properties of different porous materials are calculated through a homogenization method, macroscopic layout of different microstructures is optimized through an ordered SIMP method, designable connection domains with the same topological description are arranged in boundary areas of the microstructures of the different porous materials, and consistency sensitivity calculation of a space-time discrete power system is realized based on a first discrete-last differential adjoint sensitivity analysis method. The multi-microstructure multi-scale parallel topology optimization method under the action of the time domain dynamic load comprises the following steps: defining an initial design domain and optimization algorithm parameters, calculating equivalent elastic matrixes and mass densities of different porous materials on a microscopic scale by utilizing a microstructure finite element equation and a homogenization method, constructing a material model, and identifying structural domains and interpolation stiffness matrixes of different porous materials; performing macroscopic-scale transient finite element analysis by a HHT-alpha method, calculating an objective function and a constraint function, and obtaining a sensitivity value of a macro/micro design variable based on a discrete-first differential accompanying sensitivity analysis method; and (3) realizing parallel iterative updating of the macroscopic structure and the microscopic structure by using an MMA method, and obtaining the multi-scale optimized topological structure of the required shape. The method is used for solving the problem of time domain dynamic load optimization in the field of engineering application.

Description

Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain
Technical Field
The invention relates to the field of engineering structure optimization design, in particular to a multi-microstructure multi-scale parallel topology optimization method and system under the action of time domain dynamic load.
Background
The porous material has large porosity, light weight and good multiple physical properties, such as shock absorption, energy absorption, heat insulation, noise reduction and the like, and is widely applied to the fields of aerospace and the like. The multi-scale structure optimization design fully exploits the application potential of materials from macroscopic and microscopic simultaneously, is an important way for developing high-performance, lightweight and multifunctional advanced structures, and is widely focused and researched by students at home and abroad. However, the current research results are all based on a multi-scale uniform structure of a single unit cell, and neglect the influence of the differences of the bearing states of different material points of the structure on the configuration of the unit cell, so that the bearing capacity of the unit cell is relatively poor.
The lattice structure can bear different load, and the material characteristics of various unit cell structures or single crystal cell lattice structures with different gradient distribution can be fully utilized, so that the multi-physical properties of the composite structure can be improved. Numerous scholars consider constructing a series of multi-cell structural designs starting from the geometric continuity of adjacent microstructures, but neglecting connectivity of different cell structures in different macroscopic sub-domains, thereby restricting manufacturability of multi-scale optimization. Research work on multi-cell and gradient lattice structures has focused on optimization fields such as linear statics, frequency domain characteristics and frequency domain dynamic response, while the problem of time domain dynamic topology optimization involves complex transient sensitivity analysis, resulting in relatively less related research work.
Typically, sensitivity analysis involves only differential operations on design variables, but the time domain response sensitivity problem also involves discretization of the time domain. Thus, differentiation and discrete sequencing may have an impact on the time response sensitivity results. The sensitivity analysis strategy of firstly differentiation and then dispersion is to apply differential operation to a dynamic finite element model with continuous time domain, and then obtain the structural response field of each discrete moment through a numerical algorithm; the sensitivity analysis strategy of the discrete-first differential method is to construct the accompanying equation on a time-discrete dynamic finite element model. Although a differential-then-discrete sensitivity analysis strategy is easier to implement, such sensitivity analysis strategy can lead to consistency errors that can result in inaccuracy in the optimization results, while employing a differential-then-discrete concomitant variable approach can greatly reduce consistency errors.
The prior art with the document number of CN116150834B discloses a parallel topological optimization method (patent number of 202211224295.1) for the problem of dynamic stiffness in the time domain of a double-scale hierarchical structure, which is mainly multi-scale optimization of a single material and focuses on sensitivity analysis, and connectivity factors of different unit cell structures in different macroscopic sub-domains are not considered.
Disclosure of Invention
The technical problems to be solved by the invention are as follows:
aiming at the defects and shortcomings of the problems, the invention provides a multi-microstructure multi-scale parallel topology optimization method and a multi-scale parallel topology optimization system under the action of time domain dynamic load.
The invention firstly introduces a double Helmholtz smooth-blocking projection method to identify macroscopic distribution areas of different unit cell microstructures, calculates equivalent elastic matrixes and equivalent density parameters of the double-scale structure by using a numerical homogenization method, establishes a dynamics multi-scale optimization model, solves transient response of the double-scale structure by adopting an HHT-alpha method, calculates sensitivity of a transient dynamics finite element model by using a discrete-differential-first-followed variable method, and finally realizes updating of the macro microstructure by using an MMA method. The invention also combines two specific numerical examples of the cantilever beam and the clamped beam to verify the effectiveness of the parallel dynamic topology optimization method of the multi-microstructure under the action of the dynamic load in the time domain.
The technical scheme adopted by the invention for solving the technical problems is as follows:
A multi-microstructure multi-scale parallel topology optimization method under the action of time domain dynamic load comprises the following steps: s1, setting a microstructure connection domain, initializing design variables of macro and micro dimensions, carrying out finite element analysis on a microstructure by adopting a homogenization method, and calculating an equivalent elastic matrix and equivalent density to form an initialization model; s2, introducing a double Helmholtz smooth-block projection scheme, and carrying out interpolation on macroscopic design variables and an equivalent elastic matrix by using an Ordered-SIMP method; s3, adopting an HHT-alpha method as a time integration scheme of a time domain dynamics model, solving transient response of the double-scale structure, and deducing sensitivity of a transient dynamics problem objective function and a constraint function based on a discrete-differential-first-associated variable method; s4, realizing parallel iterative updating of the macroscopic and microscopic structures by using a mobile asymptote method (MMA) method, and obtaining the multi-scale optimized topological structure of the required shape.
In the step S1, a microstructure connection domain is set, design variables of a macro scale and a micro scale are initialized, finite element analysis is carried out on a microstructure by adopting a homogenization method, an equivalent elastic matrix and an equivalent density are calculated, and an initialization model is formed;
on the macrostructure scale, the desired values of density of the different porous materials are arranged in ascending order Expressed by the normalization process as:
wherein: ρ max is the maximum value of the density expectations in m single cell structures; the void cells are considered as regions of very low density, and are typically arranged to avoid singularities in the stiffness matrix When the microstructures are designed, the same finite element grids are adopted to disperse each microstructure, so that the cell design variable values at the same position in the designable connection area are equal, namely:
Wherein: n c is the number of finite elements of the microstructure attachment region, Representing the cell design variables at the same location within the connection region.
In the step S2, a double Helmholtz smooth-block projection scheme is introduced, and interpolation of macroscopic design variables and an equivalent elastic matrix is carried out by using an Ordered-SIMP method;
introducing a double Helmholtz smooth-segmented projection scheme to filter the original design variable field μ into a smooth design variable field Identifying macroscopic domains of different porous materials,/>The design variables corresponding to the original design variable field and the smooth design variable field are respectively; since the PDE filter can effectively reduce the computational cost when dealing with large filter radius problems, the PDE filter is selected to replace the standard filter in the two smoothing processes; PDE filter based on Helmholtz partial differential equation realizes twice smooth processing of macroscopic density variable, and the expression is:
wherein: (representing mu e or /) ) And X (representing/>)Or/>) Is the original and smooth density variable field, r is a length scale parameter similar to the filter radius in a standard filter,/>Representing the differential operation symbol, satisfying the following expression:
Wherein: r HS is Helmholtz smooth filter radius;
the equivalent characteristic matrix of the multi-microstructure can be effectively defined by adopting an Ordered-SIMP interpolation scheme, and the optimization problem is driven to converge on a desired optimal solution by adjusting a penalty factor; order the (Mu e is a macroscopic design variable), a macroscopic scale cell elastic matrix/>Expressed as:
Wherein: a D and B D are scaling coefficients and conversion coefficients, p is a penalty factor; when (when) Then there are:
wherein: the equivalent elastic matrix of each phase is obtained by a calculation homogenization method.
In the step S3, a HHT-alpha method is adopted as a time integration scheme of a time domain dynamics model, transient response of a double-scale structure is solved, and sensitivity of a transient dynamics problem objective function and a constraint function is deduced based on a discrete-differential-first-followed variable method;
The dynamic compliance minimization problem is defined as:
find:
min:
s.t.:
Wherein: f is the target dynamic compliance, M and K are the global mass matrix and stiffness matrix, respectively, C is the damping matrix (C=α rM+βrK,αr、βr is the Rayleigh damping factor), And u t represent the acceleration, velocity and displacement response fields, respectively,/>, of the structure at the t-th time step under the action of the dynamic load f t A time step which is a termination time; g MA is a macro-scale volume constraint,/>Volume constraints for the ith microstructure; v MA、Vi MI is the volume usage of macro-micro material, respectively,/>Maximum volume usage allowed for macrostructure,/>Expected values for material density for ascending order; mu e is a macroscopic design variable,/>Is a microscopic design variable;
the HHT- α method modifies the finite element equation in semi-discrete form in the optimization model to:
the update format of the displacement and velocity fields is as follows by the Newmark-beta finite difference relation:
wherein: Δt represents a time relative change value, β= (1+α) 2/4, γ= (1+2α)/2 is an algorithm parameter, and reasonable selection of the numerical damping parameter α ensures that the algorithm has at least second-order accuracy and unconditional stability, and Δt represents a time increment;
the sensitivity of the objective function to macroscopic design variables is expressed by the chain derivative rule as:
for the design variables of the microstructure junction region, the sensitivity of the objective function is expressed as:
In step S4, a mobile asymptote method (MMA) is used to implement parallel iterative updating of macroscopic and microscopic structures, resulting in a multi-scale optimized topology of the desired shape, the mobile asymptote method (MMA) being a prior art category.
The system is provided with a program module corresponding to the steps of the technical scheme, and the steps in the multi-microstructure multi-scale parallel topology optimization method under the action of the time domain dynamic load are executed during operation. A computer readable storage medium storing a computer program configured to implement the steps of the above-described method for multi-microstructure multi-scale parallel topology optimization under time domain dynamic loading when called by a processor.
The invention has the following beneficial technical effects:
The structural design of the multicrystal cell specifically considers the connectivity of different cell structures in different macroscopic domains so as to promote the manufacturability of multi-scale optimization.
According to the method, a double Helmholtz smooth-blocking projection scheme is introduced, macroscopic structural domains of different porous materials are identified, macroscopic equivalent mechanical properties of the porous materials are calculated through a homogenization method, macroscopic layout of different microstructures is optimized through an ordered SIMP method, designable connection domains with the same topological description are arranged in boundary areas of the microstructures of the different porous materials, and consistency sensitivity calculation of a space-time discrete power system is realized based on a first discrete-then differential adjoint sensitivity analysis method. The multi-microstructure multi-scale parallel topology optimization method under the action of the time domain dynamic load comprises the following steps: defining an initial design domain and optimization algorithm parameters, calculating equivalent elastic matrixes and mass densities of different porous materials on a microscopic scale by utilizing a microstructure finite element equation and a homogenization method, constructing a material model, and identifying structural domains and interpolation stiffness matrixes of different porous materials; performing macroscopic-scale transient finite element analysis by a HHT-alpha method, calculating an objective function and a constraint function, and obtaining a sensitivity value of a macro/micro design variable based on a discrete-first differential accompanying sensitivity analysis method; and (3) realizing parallel iterative updating of the macroscopic structure and the microscopic structure by using an MMA method, and obtaining the multi-scale optimized topological structure of the required shape.
The method is used for solving the problem of optimizing the dynamic load of the time domain in the field of engineering application.
Drawings
Fig. 1 is a flowchart of a multi-microstructure multi-scale parallel topology optimization method under the action of a time domain dynamic load provided by an embodiment of the invention.
FIG. 2 is a process diagram of a segmented projection method employed by the present invention to effect identification of macroscopic regions of different microstructure blocks.
FIG. 3 is a segmented projection scheme of different porous materials employed in the present invention.
FIG. 4 is a diagram of a programmable attachment region approach taken by the present invention to ensure connectivity between microstructures of different porous materials.
FIG. 5 is a schematic diagram of the design domain, load, and boundary conditions of the cantilever beam provided by the present invention.
Fig. 6 is a schematic diagram of sinusoidal loading of cantilever Liang Banbo provided by the present invention.
FIG. 7 shows sensitivity calculations for cantilever Liang Hongguan unit (a) and cell 1 (b), cell 2 (c), and cell 3 (d) provided by the present invention using first discrete-then differential.
Fig. 8 shows a cantilever Liang Mubiao function and constraint function convergence procedure (a) and an optimization topology (b) provided by the present invention.
FIG. 9 shows the optimized results of the cantilever beam provided by the invention at load time of 0.03s (a) and 0.01s (b).
FIG. 10 is a graph of the vertical displacement versus structural transient compliance history for cantilever load points provided by the present invention at 0.05s (a), 0.03s (b) and 0.01s (c) of time of action.
FIG. 11 is a schematic diagram of the design domain, load, and boundary conditions of the clamped beam provided by the present invention.
Fig. 12 is a schematic diagram of half-wave cosine dynamic load of a clamped beam provided by the invention.
FIG. 13 shows the optimized results of the clamped beams provided by the invention at load application times of 0.5s (a), 0.1s (b) and 0.03s (c).
FIG. 14 shows the vertical displacement and structure transient compliance history of the clamped beam load point provided by the present invention at 0.5s (a), 0.1s (b) and 0.03s (c) of time of action.
FIG. 15 shows the macro and micro optimization design results of the single cell clamped beams provided by the present invention.
FIG. 16 shows the macro and micro optimization design results of 5 unit cell clamped beams provided by the invention.
Detailed Description
1-16, The implementation of the multi-microstructure multi-scale parallel topology optimization method under the action of the time domain dynamic load is described as follows:
as shown in fig. 1, the flow of the method of the present invention is as follows:
S1, setting a microstructure connection domain, initializing design variables of macro and micro dimensions, carrying out finite element analysis on a microstructure by adopting a homogenization method, and calculating an equivalent elastic matrix and equivalent density to form an initialization model.
S2, introducing a double Helmholtz smooth-block projection scheme, and carrying out interpolation on macroscopic design variables and an equivalent elastic matrix by using an Ordered-SIMP method;
S3, adopting an HHT-alpha method as a time integration scheme of a time domain dynamics model, solving transient response of the double-scale structure, and deducing sensitivity of a transient dynamics problem objective function and a constraint function based on a discrete-differential-first-adjoint variable method.
S4, realizing parallel iterative updating of the macroscopic and microscopic structures by using a mobile asymptote method (MMA) method, and obtaining the multi-scale optimized topological structure of the required shape.
Further, in step S1, in order to ensure connectivity between microstructures of different porous materials, a designable connection region Ω c with a designated width is disposed near the boundary of the microstructures, in which each microstructure is used to share the same micro design variable and filter radius, and mutually independent micro design variables Ω q and Ω w are used in the inner region of the microstructure. When the microstructures are designed, the same finite element grids are adopted to disperse each microstructure, so that the cell design variable values at the same position in the designable connection area are equal, namely:
Wherein: n c is the number of finite elements of the microstructure attachment region, Respectively representing cell design variable values at the same location. Different microstructures have the same topological description on their boundaries, enhancing manufacturability of the optimum design results. In addition, the method can effectively ensure good connectivity of the microstructure designs of porous materials with different configurations under the condition that optimization constraint conditions are not applied, as shown in fig. 4.
On the microscopic scale, introducing m kinds of porous materials with microstructure characterization different, and obtaining the equivalent elastic matrix of each phase by a calculation homogenization methodOn a macroscopic scale, global volume constraints are utilized to limit the volumetric usage of each material phase. Subsequently, the elastic matrix/>, of the macro-cell is calculated using an interpolation strategyThe dynamic compliance of the structure is thus given based on dynamic finite element analysis.
For m porous materials, define a density design variableThe i-th microstructure topology was characterized (i=1, 2,..m; j=1, 2,..n MI), where n MI represents the finite number of units for the discrete porous material microstructure unit cell. The macro-scale density design variable μ e(e=1,2,...,nMA), where n MA represents the finite number of cells for a discrete macrostructure. The dynamic compliance minimization problem is defined as:
find:
min:
s.t.:
Wherein: f is the target dynamic compliance, M and K are the global mass matrix and stiffness matrix, respectively, C is the damping matrix (C=α rM+βrK,αr、βr is the Rayleigh damping factor), And u t represent the acceleration, velocity and displacement response fields, respectively,/>, of the structure at the t-th time step under the action of the dynamic load f t Is the time step of the termination moment. G MA is a macro-scale volume constraint,/>Is the volume constraint of the ith microstructure. /(I)Maximum volume usage allowed for macrostructure,/>Is the desired value of the density of the porous material. The macro/micro material volume usage V MA、Vi MI is obtained by:
wherein: And/> The unit volumes at macro/micro scale, respectively.
Given that the unit cell size is much smaller than the macrostructure size, the scale separable assumption is satisfied; in addition, the unit cells are uniformly distributed on a macroscopic scale, satisfying periodic boundary conditions. Thus, based on the calculation homogenization method, the equivalent elastic matrixExpressed as:
wherein: representing the stress-strain matrix of the porous material i in the planar stress state,/> And/>Respectively, an elastic model and poisson ratio, epsilon 0 is a unit test strain field which is independent of linearity, and epsilon (lambda i) and lambda i are strain fields and displacement fields under different unit strain.
Lambda i is divided into 4 parts: the specified displacement lambda i1 and the unknown displacement lambda i2 of the internal node meet the displacement of the periodic boundary condition(/>To the differential of the displacement of the microstructure across the edge). Applying periodic boundary conditions, the equilibrium equation for the unit cell is:
Wherein: f 1 is a constraint counter force, and the stiffness matrix K MI is a symmetric matrix according to the periodic boundary conditions F 2 =0 and F 3+F4 =0.
Eliminating the first line of the balance equation, adding the third line to the fourth line, and usingThe equilibrium equation for the original unit cell is reduced to:
Similarly, using first order homogenization theory, the mass density of unit e at the macro scale Can be expressed as:
wherein: v i MI is the volume of the microstructure of the porous material i, which is the mass density of the microstructure-based material,/> And/>The cell volume and density variations of the microstructure of the porous material i, respectively.
Further, in step S2, a double Helmholtz smoothing-chunking projection scheme is introduced to filter the original design variable field μ into a smooth design variable fieldMacroscopic domains of different porous materials are identified. Interpolation of the macroscopic design variables and the equivalent elastic matrix is performed using the Ordered-SIMP method. The specific process of the double Helmholtz smooth-block projection scheme is shown in figures 2 and 3;
on the macrostructure scale, the desired values of density of the different porous materials are arranged in ascending order Expressed by the normalization process as:
Wherein: ρ max is the maximum value of the density expectations in m unit cell structures, Is a different porous material expected value for the normalization process. The void cells are considered as regions of very low density, and are typically arranged to avoid singularities in the stiffness matrixWill/>Is introduced into a standard smooth filter and is expanded into a segmented projection filter made of different porous materials. For the followingThen there are:
wherein: (representative/> Or/>) Is the unit density variable to be filtered,/>(Representative/>Or/>) Is a block projection filtered unit density variable, x A represents an intermediate interpolation density variable. Control threshold parameter η pp and sharpness parameter β pp Extrapolated projection to/>Or/>
Since PDE filters can effectively reduce computational costs when dealing with large filter radius problems, PDE filters are chosen to replace standard filters in the two smoothing process. PDE filter based on Helmholtz partial differential equation realizes twice smooth processing of macroscopic density variable, and the expression is:
wherein: (representing mu e or /) ) And X (representing/>)Or/>) Is the original and smooth density variable field, r is a length scale parameter similar to the filter radius in a standard filter,/>Representing the differential operation symbol, satisfying the following expression:
wherein: r HS is Helmholtz smooth filter radius. The finite element method is adopted to calculate Helmholtz partial differential equation, and the following steps are:
Kf=∑eKf,e
T=∑eTe
wherein: k f is the overall smooth matrix, T is the transformation matrix, N e is the cell-shaped function array, K f,e is the cell smooth matrix, and T e is the cell transformation matrix.
The adoption of the ordered SIMP interpolation scheme can effectively define the equivalent characteristic matrix of the multi-microstructure, and the optimization problem is driven to converge on the expected optimal solution by adjusting the penalty factors.
Order the(Mu e is a macroscopic design variable), a macroscopic scale cell elastic matrix/>Expressed as:
Wherein: a D and B D are scaling and conversion coefficients, p is a penalty factor. When (when) Then there are:
In addition, for a linear dynamics system, the stiffness matrix k e and the mass matrix m e of the macro-unit are expressed as:
Further, in step S3, the HHT- α method modifies the finite element equation in semi-discrete form in the optimization model to:
the update format of the displacement and velocity fields is as follows by the Newmark-beta finite difference relation:
Wherein: Δt represents a time relative change value, beta= (1+alpha) 2/4, gamma= (1+2 alpha)/2 is an algorithm parameter, and reasonable selection of the numerical damping parameter alpha ensures that the algorithm has at least second-order precision and unconditional stability.
Substituting the displacement and velocity field update formats into the finite element equation to obtain the residual form of the momentum control equation:
Wherein: m 1、M0、C0 represents intermediate parameters of the mass matrix and the damping matrix, respectively.
For the HHT-alpha residual equation, the expression of the initial time is:
and according to the response fields such as acceleration, speed and displacement at the previous moment, iteratively updating to obtain the response fields such as acceleration, speed and displacement at the current moment.
In order to avoid consistency errors caused by a first-differentiation-last-differentiation method, a first-differentiation-last-differentiation method is adopted to complete sensitivity analysis of the linear power system. Assuming that the load and the initial conditions are independent of the design variables, x is used for replacing macro-micro design variables mu e andThe sensitivity of the objective function to the design variable is expressed as:
The update format of the displacement and velocity field obtained by the calculation of the Newmark-beta finite difference relation is converted into the residual format, and then the following steps are:
Wherein: p t、Qt is the residual format of the displacement, velocity field, respectively.
The sensitivity of the objective function to the design variables is rewritten as follows, introducing the companion variable λ t、μt、ζt:
Since P t and Q t do not contain design variables, there are The initial conditions are independent of the design variables, then there is/>Elimination of inclusion/>And/>The items of (1) are:
finally, the sensitivity of the objective function to the design variables is reduced to:
Substituting the HHT-alpha residual equation, the initial condition and the Newmark-beta finite difference relation into the accompanying equation to obtain solutions of all accompanying variables:
ζt-1=C0λt+Δtμtt
the sensitivity of the objective function to the macroscopic design variable μ e by the chain derivative rule is expressed as:
Objective function for And/>The sensitivity of (c) is expressed as:
Since the objective function is not apparent And/>Then there are:
in addition, the residual error of the momentum control equation and the HHT-alpha residual error equation of the initial moment are respectively matched And/>The partial derivative is calculated by:
Macroscopic scale cell elastic matrix And macroscopic volume constraint G MA sensitivity with respect to μ e are respectively:
To facilitate programming, the design variables for each microstructure are re-indexed: connection region Change to/>Internal region/>By means of the chain-derivative principle, the objective function is designed for internal microscopic design variables/>The sensitivity of (c) is expressed as:
The sensitivity of the objective function to each equivalent elastic matrix D HMI and equivalent density ρ HMI resulting from the homogenization method is expressed as:
In addition, the residual error of the momentum control equation and the HHT-alpha residual equation at the initial moment are respectively matched with D HMI and D HMI The partial derivative is calculated by:
/>
for the design variables of the microstructure junction region, the sensitivity of the objective function is expressed as:
Finally, microscopic volume constraints For/>The sensitivity of (2) is:
Further, in step S4, a mobile asymptote method (MMA) method is selected, so that the method has better adaptability, and can solve a wider and more complex optimization problem. The mobile asymptote method (MMA) resolves an implicit problem into a plurality of explicit convex approximation sub-problems. Thereby conveniently realizing parallel iterative updating of macroscopic and microscopic design variables, obtaining a multi-scale optimized topological structure of a required shape, and determining an optimal solution of an optimization problem according to a convergence criterion:
Wherein: ζ is the allowable convergence error limit.
The present invention will be described in further detail with reference to the drawings and examples, in order to make the objects, technical solutions and advantages of the present invention more apparent. It should be understood that the specific embodiments described herein are for purposes of illustration only and are not intended to limit the scope of the invention. In addition, technical features of the embodiments of the present invention described below may be combined with each other as long as they do not collide with each other.
Referring to fig. 5-16, the present invention verifies the effectiveness of the optimization algorithm under different dynamic load input conditions for both cantilever and clamped beams models.
The invention provides a multi-microstructure multi-scale parallel topology optimization method under the action of time domain dynamic load, which comprises the following steps:
S1, setting a microstructure connection domain, initializing design variables of macro and micro dimensions, carrying out finite element analysis on a microstructure by adopting a homogenization method, and calculating an equivalent elastic matrix and equivalent density to form an initialization model.
Referring to fig. 5 and 6, the cantilever beam adopts the dimension parameters of length l=8m, width h=4m and thickness h=0.01 m, and selects steel base material and elastic modulusPoisson's ratio/>Mass density/>The macro design domain and the micro design domain are respectively discretized into 20000 and 2500 bilinear rectangular units. The Rayleigh damping parameter of the base material is alpha r=10s-1r=10-5 s, and the number of response sampling points/>Macroscopic volume constraint parameter/>Microcosmic volume constraint parameter/>Set to [10 -9, 0.5,0.5,0.5].
Referring to fig. 11 and 12, the clamped beam adopts dimension parameters of length l=12m, width h=2m, and thickness h=0.01 m, and selects steel base material, elastic modulusPoisson's ratio/>Mass density/>The macro design domain and the micro design domain are respectively discretized into 20000 and 2500 bilinear rectangular units. The Rayleigh damping parameter of the base material is alpha r=10s-1r=10-5 s, and the number of response sampling points/>Macroscopic volume constraint parameter/>Microcosmic volume constraint parameter/>Set to [10 -9, 0.4,0.4,0.4].
In order to ensure connectivity between the microstructures of different porous materials, a programmable connection region omega c with a specified width is arranged near the boundary of the microstructures, the microstructures are adopted to share the same micro design variable and filter radius in the region, and mutually independent micro design variables omega q and omega w are adopted in the inner region of the microstructures. When the microstructures are designed, the same finite element grids are adopted to disperse each microstructure, so that the cell design variable values at the same position in the designable connection area are equal, namely:
Wherein: n c is the number of finite elements of the microstructure attachment region, Respectively representing cell design variable values at the same location. Different microstructures have the same topological description on their boundaries, enhancing manufacturability of the optimum design results. In addition, the method can effectively ensure good connectivity of the microstructure design of porous materials with different configurations under the condition that optimization constraint conditions are not applied.
On the microscopic scale, introducing m kinds of porous materials with microstructure characterization different, and obtaining the equivalent elastic matrix of each phase by a calculation homogenization methodOn a macroscopic scale, global volume constraints are utilized to limit the volumetric usage of each material phase. Subsequently, the elastic matrix/>, of the macro-cell is calculated using an interpolation strategyThe dynamic compliance of the structure is thus given based on dynamic finite element analysis.
For m porous materials, define a density design variableThe i-th microstructure topology was characterized (i=1, 2,..m; j=1, 2,..n MI), where n MI represents the finite number of units for the discrete porous material microstructure unit cell. The macro-scale density design variable μ e(e=1,2,...,nMA), where n MA represents the finite number of cells for a discrete macrostructure. The dynamic compliance minimization problem is defined as:
find:
min:
s.t.:
Wherein: f is the target dynamic compliance, M and K are the global mass matrix and stiffness matrix, respectively, C is the damping matrix (C=α rM+βrK,αr、βr is the Rayleigh damping factor), And u t represent the acceleration, velocity and displacement response fields, respectively,/>, of the structure at the t-th time step under the action of the dynamic load f t Is the time step of the termination moment. G MA is a macro-scale volume constraint,/>Is the volume constraint of the ith microstructure. /(I)Maximum volume usage allowed for macrostructure,/>Is the desired value of the density of the porous material. The macro/micro material volume usage V MA、Vi MI is obtained by:
/>
wherein: And/> The unit volumes at macro/micro scale, respectively.
Given that the unit cell size is much smaller than the macrostructure size, the scale separable assumption is satisfied; in addition, the unit cells are uniformly distributed on a macroscopic scale, satisfying periodic boundary conditions. Thus, based on the calculation homogenization method, the equivalent elastic matrixExpressed as:
wherein: representing the stress-strain matrix of the porous material i in the planar stress state,/> And/>Respectively, an elastic model and poisson ratio, epsilon 0 is a unit test strain field which is independent of linearity, and epsilon (lambda i) and lambda i are strain fields and displacement fields under different unit strain.
Lambda i is divided into 4 parts: the specified displacement lambda i1 and the unknown displacement lambda i2 of the internal node meet the displacement of the periodic boundary condition(/>To the differential of the displacement of the microstructure across the edge). Applying periodic boundary conditions, the equilibrium equation for the unit cell is:
Wherein: f 1 is a constraint counter force, and the stiffness matrix K MI is a symmetric matrix according to the periodic boundary conditions F 2 =0 and F 3+F4 =0.
Eliminating the first line of the balance equation, adding the third line to the fourth line, and usingThe equilibrium equation for the original unit cell is reduced to:
Similarly, using first order homogenization theory, the mass density of unit e at the macro scale Can be expressed as:
wherein: v i MI is the volume of the microstructure of the porous material i, which is the mass density of the microstructure-based material,/> And/>The cell volume and density variations of the microstructure of the porous material i, respectively.
S2, introducing a double Helmholtz smooth-block projection scheme, and filtering an original design variable field mu into a smooth design variable fieldMacroscopic domains of different porous materials are identified. Interpolation of the macroscopic design variables and the equivalent elastic matrix is performed using the Ordered-SIMP method.
On the macrostructure scale, the desired values of density of the different porous materials are arranged in ascending orderExpressed by the normalization process as:
Wherein: ρ max is the maximum value of the density expectations in m unit cell structures, Is a different porous material expected value for the normalization process. The void cells are considered as regions of very low density, and are typically arranged to avoid singularities in the stiffness matrixWill/>Is introduced into a standard smooth filter and is expanded into a segmented projection filter made of different porous materials. For the followingThen there are:
wherein: (representative/> Or/>) Is the unit density variable to be filtered,/>(Representative/>Or/>) Is a block projection filtered unit density variable, x A represents an intermediate interpolation density variable. Control threshold parameter η pp and sharpness parameter β pp Extrapolated projection to/>Or/>
Since PDE filters can effectively reduce computational costs when dealing with large filter radius problems, PDE filters are chosen to replace standard filters in the two smoothing process. PDE filter based on Helmholtz partial differential equation realizes twice smooth processing of macroscopic density variable, and the expression is:
wherein: (representing mu e or /) ) And X (representing/>)Or/>) Is the original and smooth density variable field, r is a length scale parameter similar to the filter radius in a standard filter,/>Representing the differential operation symbol, satisfying the following expression:
wherein: r HS is Helmholtz smooth filter radius. The finite element method is adopted to calculate Helmholtz partial differential equation, and the following steps are:
Kf=∑eKf,e
T=∑eTe
wherein: k f is the overall smooth matrix, T is the transformation matrix, N e is the cell-shaped function array, K f,e is the cell smooth matrix, and T e is the cell transformation matrix.
The adoption of the ordered SIMP interpolation scheme can effectively define the equivalent characteristic matrix of the multi-microstructure, and the optimization problem is driven to converge on the expected optimal solution by adjusting the penalty factors.
Order the(Mu e is a macroscopic design variable), a macroscopic scale cell elastic matrix/>Expressed as:
Wherein: a D and B D are scaling and conversion coefficients, p is a penalty factor. When (when) Then there are:
In addition, for a linear dynamics system, the stiffness matrix k e and the mass matrix m e of the macro-unit are expressed as:
S3, adopting an HHT-alpha method as a time integration scheme of a time domain dynamics model, solving transient response of the double-scale structure, and deducing sensitivity of a transient dynamics problem objective function and a constraint function based on a discrete-differential-first-adjoint variable method.
The HHT- α method modifies the finite element equation in semi-discrete form in the optimization model to:
the update format of the displacement and velocity fields is as follows by the Newmark-beta finite difference relation:
Wherein: Δt represents a time relative change value, beta= (1+alpha) 2/4, gamma= (1+2 alpha)/2 is an algorithm parameter, and reasonable selection of the numerical damping parameter alpha ensures that the algorithm has at least second-order precision and unconditional stability.
Substituting the displacement and velocity field update formats into the finite element equation to obtain the residual form of the momentum control equation:
Wherein: m 1、M0、C0 represents intermediate parameters of the mass matrix and the damping matrix, respectively.
For the HHT-alpha residual equation, the expression of the initial time is:
and according to the response fields such as acceleration, speed and displacement at the previous moment, iteratively updating to obtain the response fields such as acceleration, speed and displacement at the current moment.
In order to avoid consistency errors caused by a first-differentiation-last-differentiation method, a first-differentiation-last-differentiation method is adopted to complete sensitivity analysis of the linear power system. Assuming that the load and the initial conditions are independent of the design variables, x is used for replacing macro-micro design variables mu e andThe sensitivity of the objective function to the design variable is expressed as: /(I)
The update format of the displacement and velocity field obtained by the calculation of the Newmark-beta finite difference relation is converted into the residual format, and then the following steps are:
Wherein: p t、Qt is the residual format of the displacement, velocity field, respectively.
The sensitivity of the objective function to the design variables is rewritten as follows, introducing the companion variable λ t、μt、ζt:
Since P t and Q t do not contain design variables, there are The initial conditions are independent of the design variables, then there is/>Elimination of inclusion/>And/>The items of (1) are:
finally, the sensitivity of the objective function to the design variables is reduced to:
Substituting the HHT-alpha residual equation, the initial condition and the Newmark-beta finite difference relation into the accompanying equation to obtain solutions of all accompanying variables:
ζt-1=C0λt+Δtμtt
the sensitivity of the objective function to the macroscopic design variable μ e by the chain derivative rule is expressed as:
Objective function for And/>The sensitivity of (c) is expressed as: /(I)
Since the objective function is not apparentAnd/>Then there are:
in addition, the residual error of the momentum control equation and the HHT-alpha residual error equation of the initial moment are respectively matched And/>The partial derivative is calculated by:
Macroscopic scale cell elastic matrix And macroscopic volume constraint G MA sensitivity with respect to μ e are respectively:
To facilitate programming, the design variables for each microstructure are re-indexed: connection region Change to/>Internal region/>By means of the chain-derivative principle, the objective function is designed for internal microscopic design variables/>The sensitivity of (c) is expressed as:
The sensitivity of the objective function to each equivalent elastic matrix D HMI and equivalent density ρ HMI resulting from the homogenization method is expressed as:
/>
In addition, the residual error of the momentum control equation and the HHT-alpha residual equation at the initial moment are respectively matched with D HMI and D HMI The partial derivative is calculated by:
for the design variables of the microstructure junction region, the sensitivity of the objective function is expressed as:
Finally, microscopic volume constraints For/>The sensitivity of (2) is:
S4, realizing parallel iterative updating of the macroscopic and microscopic structures by using a mobile asymptote method (MMA) method, and obtaining the multi-scale optimized topological structure of the required shape.
The mobile asymptote method (MMA) method is selected, so that the method has better adaptability, and can solve a wider and more complex optimization problem. The mobile asymptote method (MMA) resolves an implicit problem into a plurality of explicit convex approximation sub-problems. Thereby conveniently realizing parallel iterative updating of macroscopic and microscopic design variables, obtaining a multi-scale optimized topological structure of a required shape, and determining an optimal solution of an optimization problem according to a convergence criterion:
Wherein: ζ is the allowable convergence error limit.
The optimization result obtained by the cantilever Liang Shili adopted by the invention shows that: the sensitivity calculation value of the cantilever beam macroscopic unit and the 3-class unit adopts the relative error result of the sensitivity calculation value of the discrete-differential strategy relative to the finite difference method theoretical value, and the effectiveness of the accompanying sensitivity analysis method of the transient dynamics optimization problem of the multicrystal cell double-scale structure is verified; different porous materials have distinct topologies, however, the optimized unit cell structure can be effectively connected in the horizontal and vertical directions, effectively enhancing the manufacturability of the dual-scale structural topology optimization design.
The optimization result obtained by the clamped beam example adopted by the invention shows that: when a load is rapidly applied to the structure, the optimizer places more reinforcing ribs in the central region of the structure to counter the progressively increasing dynamic effects. The amplitude attenuation of the vertical displacement and the transient compliance of the structure is obvious due to the dissipation effect of damping, and the reduction of the load acting time leads to the weakening of the effect of the damping effect on the dynamic response of the structure. With the increase of the variety of porous materials, the design space of the structure-material integration is more abundant, so that the numerical value of the optimal dynamic flexibility is obviously reduced, and the structural dynamic stiffness is further improved.
The results of two examples of the present invention show that: the method combines the identification algorithm of different porous material blocks and the method of a designable connection area for ensuring the connectivity of the microstructure, so that the integrated topological optimization of the structure and the material containing various porous materials under the transient excitation effect can be realized; the accompanying sensitivity analysis method of the transient dynamics optimization problem of the multicrystal cell double-scale structure based on the discrete-differential strategy effectively ensures the reliable solution of the multi-scale parallel optimization problem; in addition, with the increase of the types of porous materials, the design space of the structural optimization problem of the composite material is prolonged, and the target performance is further improved.

Claims (3)

1. A multi-microstructure multi-scale parallel topology optimization method under the action of time domain dynamic load is characterized by comprising the following steps of: the method comprises the following steps:
s1, setting a microstructure connection domain, initializing design variables of macro and micro dimensions, carrying out finite element analysis on a microstructure by adopting a homogenization method, and calculating an equivalent elastic matrix and equivalent density to form an initialization model;
s2, introducing a double Helmholtz smooth-block projection scheme, and carrying out interpolation on macroscopic design variables and an equivalent elastic matrix by using an Ordered-SIMP method;
S3, adopting an HHT-alpha method as a time integration scheme of a time domain dynamics model, solving transient response of the double-scale structure, and deducing sensitivity of a transient dynamics problem objective function and a constraint function based on a discrete-differential-first-associated variable method;
S4, implementing parallel iterative updating of the macroscopic structure and the microscopic structure by using a moving asymptote method to obtain a multi-scale optimized topological structure with a required shape;
the implementation process of the step S1 is as follows:
On the macrostructure scale, the density expectations of the different porous materials are arranged in ascending order as ρ i AO, expressed by a normalization process:
wherein: ρ max is the maximum value of the density expectations in m single cell structures; the void cells are considered as regions of very low density, and are typically arranged to avoid singularities in the stiffness matrix
When the microstructures are designed, the same finite element grids are adopted to disperse each microstructure, so that the cell design variable values at the same position in the designable connection area are equal, namely:
Wherein: n c is the number of finite elements of the microstructure attachment region, Respectively representing unit design variables at the same position in the connection area;
The implementation process of the step S2 is as follows:
introducing a double Helmholtz smooth-segmented projection scheme to filter the original design variable field μ into a smooth design variable field Identifying macroscopic domains of different porous materials,/>The design variables corresponding to the original design variable field and the smooth design variable field are respectively; since the PDE filter can effectively reduce the computational cost when dealing with large filter radius problems, the PDE filter is selected to replace the standard filter in the two smoothing processes;
PDE filter based on Helmholtz partial differential equation realizes twice smooth processing of macroscopic density variable, and the expression is:
wherein: (representing mu e or /) ) And X (representing/>)Or/>) Is the original and smooth density variable field, r is a length scale parameter similar to the filter radius in a standard filter,/>Representing the differential operation symbol, satisfying the following expression:
Wherein: r HS is Helmholtz smooth filter radius;
the equivalent characteristic matrix of the multi-microstructure can be effectively defined by adopting an Ordered-SIMP interpolation scheme, and the optimization problem is driven to converge on a desired optimal solution by adjusting a penalty factor; order the (Mu e is a macroscopic design variable), a macroscopic scale cell elastic matrix/>Expressed as:
Wherein: a D and B D are scaling coefficients and conversion coefficients, p is a penalty factor; when (when) Then there are:
wherein: Obtaining an equivalent elastic matrix of each phase through a calculation homogenization method;
In the step S3, a HHT-alpha method is adopted as a time integration scheme of a time domain dynamics model, transient response of a double-scale structure is solved, and sensitivity of a transient dynamics problem objective function and a constraint function is deduced based on a discrete-differential-first-followed variable method;
The dynamic compliance minimization problem is defined as:
wherein: f is the target dynamic compliance, M and K are the global mass matrix and stiffness matrix, C is the damping matrix, C=α rM+βrK,αr、βr is the Rayleigh damping factor, And u t represent the acceleration, velocity and displacement response fields, respectively,/>, of the structure at the t-th time step under the action of the dynamic load f t A time step which is a termination time; g MA is a macro-scale volume constraint,/>Volume constraints for the ith microstructure; v MA、Vi MI is the volume usage of macro-micro material, respectively,/>Maximum volume usage allowed for macrostructure,/>Expected values for material density for ascending order; mu e is a macroscopic design variable,/>Is a microscopic design variable;
the HHT- α method modifies the finite element equation in semi-discrete form in the optimization model to:
the update format of the displacement and velocity fields is as follows by the Newmark-beta finite difference relation:
wherein: Δt represents a time relative change value, β= (1+α) 2/4, γ= (1+2α)/2 is an algorithm parameter, and reasonable selection of the numerical damping parameter α ensures that the algorithm has at least second-order accuracy and unconditional stability, and Δt represents a time increment;
the sensitivity of the objective function to macroscopic design variables is expressed by the chain derivative rule as:
for the design variables of the microstructure junction region, the sensitivity of the objective function is expressed as:
2. A multi-microstructure multi-scale parallel topology optimization system under the action of time domain dynamic load is characterized in that: the system has program modules corresponding to the steps of any of the preceding claims 1, and the steps of the multi-microstructure multi-scale parallel topology optimization method under the action of a time domain dynamic load are executed in the running process.
3. A computer-readable storage medium, characterized by: the computer readable storage medium stores a computer program configured to implement the steps of the multi-microstructure multi-scale parallel topology optimization method under a time domain dynamic load of any of claims 1 when invoked by a processor.
CN202311492300.1A 2023-11-09 Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain Active CN117373579B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202311492300.1A CN117373579B (en) 2023-11-09 Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202311492300.1A CN117373579B (en) 2023-11-09 Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain

Publications (2)

Publication Number Publication Date
CN117373579A CN117373579A (en) 2024-01-09
CN117373579B true CN117373579B (en) 2024-05-17

Family

ID=

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112182929A (en) * 2020-09-18 2021-01-05 北京航空航天大学 Size control-considered cross-scale reliability topological optimization method for porous material
CN112417692A (en) * 2020-11-24 2021-02-26 华东交通大学 Multi-scale topological optimization design method of material structure based on load uncertainty
CN115577447A (en) * 2022-09-28 2023-01-06 北京理工大学 Unmanned aerial vehicle structure optimization method based on double-scale parallel topology optimization
CN116150834A (en) * 2022-10-08 2023-05-23 哈尔滨理工大学 Parallel topology optimization method for time domain dynamic stiffness problem of double-scale hierarchical structure

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112182929A (en) * 2020-09-18 2021-01-05 北京航空航天大学 Size control-considered cross-scale reliability topological optimization method for porous material
CN112417692A (en) * 2020-11-24 2021-02-26 华东交通大学 Multi-scale topological optimization design method of material structure based on load uncertainty
CN115577447A (en) * 2022-09-28 2023-01-06 北京理工大学 Unmanned aerial vehicle structure optimization method based on double-scale parallel topology optimization
CN116150834A (en) * 2022-10-08 2023-05-23 哈尔滨理工大学 Parallel topology optimization method for time domain dynamic stiffness problem of double-scale hierarchical structure

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
Concurrent multiscale topology optimisation towards design and additive manufacturing of biomimicking porous structures;Tian Lan etc.;Virtual and Physical Prototyping;20221231;第18卷(第1期);全文 *
Concurrent topology optimization of multiscale composite structures in Matlab;Jie Gao etc.;Structural and Multidisciplinary Optimization;20191231;全文 *
System reliability-based design optimization with interval parameters by sequential moving asymptote method;Zeng Meng etc.;Structural and Multidisciplinary Optimization;20211231;全文 *
具有极限弹性性能的多相材料微结构拓扑优化设计;郭名璧;内燃机与配件;20200229(第4期);第1节 *

Similar Documents

Publication Publication Date Title
CN110795873B (en) Cross-scale topology optimization method considering size control
CN111737835B (en) Three-period minimum curved surface-based three-dimensional porous heat dissipation structure design and optimization method
Wang et al. Concurrent design of hierarchical structures with three-dimensional parameterized lattice microstructures for additive manufacturing
Zhang et al. Multiscale concurrent topology optimization for cellular structures with multiple microstructures based on ordered SIMP interpolation
Wang et al. Concurrent two-scale topological design of multiple unit cells and structure using combined velocity field level set and density model
Amir et al. Approximate reanalysis in topology optimization
CN110008512B (en) Negative Poisson ratio lattice structure topology optimization method considering bearing characteristics
CN109657284B (en) Metamaterial-oriented equal-geometry topology optimization method
CN110941924B (en) Multi-component system integration integrated multi-scale topology optimization design method
CN113204906B (en) Multiphase material topology optimization design method and system considering structural stability
López et al. Model-based, multi-material topology optimization taking into account cost and manufacturability
CN113345536B (en) Structural topology optimization method based on extreme anisotropy lattice material
CN109670207B (en) Dynamic integrated design method for multiple porous material structures
Bai et al. An improved numerically-stable equivalent static loads (ESLs) algorithm based on energy-scaling ratio for stiffness topology optimization under crash loads
Xu et al. A B-spline multi-parameterization method for multi-material topology optimization of thermoelastic structures
CN111310377A (en) Non-probability reliability topological optimization design method for continuum structure under mixed constraint of fundamental frequency and frequency interval
Zhao et al. Stress‐constrained concurrent topology optimization of two‐scale hierarchical structures
Yan et al. Structure/material concurrent optimization of lattice materials based on extended multiscale finite element method
Wei et al. Topology optimization for design of hybrid lattice structures with multiple microstructure configurations
CN117373579B (en) Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain
CN111737908B (en) Skin-stringer structure rapid dynamic optimization design method based on dynamic load and static force equivalence
Yan et al. Multi-Material and Multiscale Topology Design Optimization of Thermoelastic Lattice Structures.
CN117373579A (en) Multi-microstructure multi-scale parallel topology optimization method and system under action of dynamic load in time domain
CN116187073A (en) Topology optimization method for anisotropic material transient heat transfer structure based on grid-free EFGM
CN112818583B (en) Equivalent dead load obtaining method, topology optimization method and system

Legal Events

Date Code Title Description
PB01 Publication
SE01 Entry into force of request for substantive examination
GR01 Patent grant