Periodic cellular structures with excellent mechanical properties widely exist in nature. A generative design and optimization method for triply periodic level surface (TPLS)-based functionally graded cellular structures is developed in this work. In the proposed method, by controlling the density distribution, the designed TPLS-based cellular structures can achieve better structural or thermal performances without increasing its weight. The proposed technique can be divided into four steps. First, the modified 3D implicit functions of the triply periodic minimal surfaces are developed to design different types of cellular structures parametrically and generate spatially graded cellular structures. Second, the numerical homogenization method is employed to calculate the elastic tensor and the thermal conductivity tensor of the cellular structures with different densities. Third, the optimal relative density distribution of the object is computed by the scaling laws of the TPLS-based cellular structures added optimization algorithm. Finally, the relative density of the numerical results of structure optimization is mapped into the modified parametric 3D implicit functions, which generates an optimum lightweight cellular structure. The optimized results are validated subjected to different design specifications. The effectiveness and robustness of the obtained structures is analyzed through finite element analysis and experiments. The results show that the functional gradient cellular structure is much stiffer and has better heat conductivity than the uniform cellular structure.

Introduction

Additive manufacturing (AM) allows designers to create and fabricate more complex shapes and topologies. However, current design for manufacturing techniques have been developed largely for conventional manufacturing processes [1]. Therefore, to take advantage of the design freedom enabled by AM processes, design for additive manufacturing (DfAM) methods need to be developed [1,2]. Representative DfAM methodologies include cellular structure design, topology optimization (TO), manufacturability analysis for AM, and other design approaches that can make use of the AM-enabled features [1]. Among them, an artificially designed cellular structure possesses many superior properties that make it a promising solution for various applications, such as improved stiffness [3], heat transfer [4,5], energy absorption [6], and so on. In addition, cellular structures exist widely in nature. For example, the trabecular structures of human bones, the skulls of birds, and the microstructures of bamboo all consist of gradient cellular structures [7]. These naturally evolved gradient structures distribute the material reasonably to the areas where performance needs to be optimized. These excellent performances include high specific stiffness, impact resistance, and high thermal conductivity. Inspired by gradient cellular structures in nature, we will focus on the design and optimization of functionally graded cellular structures (FGCS) in this work.

The traditional cellular structures are arranged by the regular cell structures or composed of random foam structures [3]. For instance, the designers choose the regular cellular structure in commercial software with a given volume fraction or generate the random cellular structure by the Voronoi algorithm. However, there are still some issues in designing cellular structures for AM. Usually, the lightweight design of the components use only uniform cellular structures for infilling, which cannot satisfy a reasonable stress distribution. Besides, although some optimization methods [8,9] control the cellular infill density through stress distribution based on finite element analysis (FEA), such optimizations cannot obtain the optimal solution. In contrast, the TO method can effectively optimize the distribution of materials in space. To further optimize the cellular infilling, the gradient cellular structures based on TO calculation design methods [10] have been rapidly developed. However, most methods simply match the stress size to the relative density of the cellular structures and the structural properties of the cellular itself are not considered when numerical calculating is done. In addition, these methods generally employ truss-based lattice structures, making it difficult to obtain smooth geometries that can be used for direct 3D printing. Therefore, the above problems need to be optimized and solved.

In this work, a generative design and optimization approach for the generation of functionally graded triply periodic level surface (TPLS)-based cellular structures is developed. It optimizes the mechanical stiffness or thermal conduction by redistributing the relative density of cellular structures in a design domain and under given boundary conditions. Notably, to generate graded cellular structures, a TPLS-based function representation (F-Rep) modeling method that can obtain geometrically continuous cellular structures is proposed. One of the distinctive benefits of adopting this modeling method is that the assigned parameters can control the shape and volume fraction distribution of a TPLS-based cellular structure easily. The reason why the TPLS-based cellular structures are chosen for this work is that they offer the following potential advantages over the strut-based lattice structures: (1) the open cell structure allows easy removal of powders or liquid resins; (2) the self-supporting helical structures that do not require any supporting material during the AM process [9,11]; (3) a high specific surface area and low relative density while maintaining a higher specific strength than regular diamond cellulars in certain circumstances; and (4) a functional gradient structure can be generated by the TPLS-based cellular structure easily by controlling the level-set parameter. Furthermore, the effective mechanical and thermal conductive properties of the TPLS-based cellular structures can be calculated via the computational homogenization method [12,13]. The obtained mechanical or thermal conductive properties are added to the topology optimization process to generate gradient cellular structures with effective physical meaning. Finally, for validation purposes, a variety of design specifications are given for optimization, and the effectiveness and robustness of the obtained structures are analyzed through FEA and further experiments.

Related Works

In the field of DfAM, many researchers focused on the optimization of geometric shapes and topology structures with specific mechanical properties. In this section, a brief review related to the design and optimization of cellular structures for AM is made.

TPLS-based cellular structures derived from the triply periodic minimal surface [14] are of special interest to many scientists. A typical characteristic is that the mathematical expression represents an equation composed of trigonometric functions [14], which have no intersecting or folded surfaces. “Triply periodic” means that the structure can be patterned in the 3D space, and “level surface” means that it can adjust the size of the separated area to control the density of the cellular structures. The main advantage of TPLS-based cellular structures is that its opening cell structure has high specific surface area, high porosity, and low relative density while maintaining excellent mechanical properties [15,16]. Furthermore, TPLS has been discovered to exist in the nature stable, such as the gyroid structure in butterfly wing scales [17] and the Schwartz Diamond structure in diamonds [18]. Previously, the stiffness and strength [15,16] and manufacturability by AM [9,11] of TPLS-based cellular structures were investigated numerically and experimentally. Besides, several researchers have employed TPLS-based cellular structures in a range of applications to enhance functionalities through their topology properties relationship. For example, they have used the gyroid structure as a supporting structure for AM [11], TPLS-based cellular structures as a bone tissue engineering scaffold [19], TPLS-based cellular structures with varying transport of heat and electricity [20], and TPLS-based composite structures with high conductivities [21], and so on.

Functionally graded cellular structure is a type of functionally graded material (FGM). FGM is mainly designed by having a variable proportion of composite materials with random or regular microstructures. FGCS typically refers to the variable-density cellular or porous structures design according to the physical and mechanical properties [2224]. Here, we concentrate on discussing FGCS design and optimization based on topology optimization. Brackett et al. [10] first proposed a solid isotropic material with penalization (SIMP)-based density mapping approach to optimize the intermediate densities to 2D lattice structures of varying volume fraction. Panesar et al. [25] used the density of topology optimization directly mapping into the density of the infilled lattice structures. Liu and Shapiro [26] proposed a method to generate a spatially varying structure of FGM with target properties from computed tomography image or topology optimization. Wu et al. [27] presented a self-supporting “rhombic cell” topology gradient infill design system for 3D printing object optimization. Tang et al. [28] proposed a bidirectional evolutionary structural optimization (BESO)-based method for truss-like lattice structure optimization. Alzahrani et al. [29] designed truss-like cellular structures using the relative density mapping (RDM) method. Daynes et al. [30] presented a combined methodology of topology optimization and size optimization. In addition, Cheng et al. [31] took the anisotropic property of cellular materials into account by the homogenization method. In addition, Wang et al. [32] used an efficient homogenization-based topology optimization method to generate a variable-density honeycombs cellular structure for natural frequency optimization. Cheng et al. [4,5] generated the graded lattice structures to optimize the heat conduction of the object based on level-set methods.

Topology optimization is a mathematical technique that optimizes the distribution of material in a specific design domain, with the goal of maximizing system performance for a given set of loads and constraints. The popular topology optimization approaches include the homogenization-based method [33], SIMP [34], BESO [35], the level-set method [36], the moving morphable component/moving morphable void (MMC/MMV) methods [3739], and various algorithms derived from them. The MMC method is used in conjunction with the homogenization method to generate infilling gradient microstructure faster [40]. The BESO algorithm is a finite element-based topology optimization approach, where inefficient material is iteratively removed from a structure while efficient material is simultaneously added to the structure [35,41]. The basic principle of the SIMP algorithm is to calculate the evolution of the relative density distribution of the elements over the design space [34,42]. For a comprehensive review of SIMP and BESO topology optimization techniques, see the recent survey articles [41,42].

Our work is inspired by these works, yet it avoids the use of truss-based lattice structures to generate structural compositions and avoids the density calculated from original topology optimization directly. On the one hand, truss-like cellular structures are limited in the functionally graded modeling method, and some topologies are limited in manufacture by AM. In comparison, it is easy to generate a continuous gradient cellular structure using the TPLS-based implicit modeling method. On the other hand, most of the traditional topology optimization methods are calculated for isotropic materials. However, few cellular structures are completely isotropic. Most of them are anisotropic. In our work, the elastic tensor and the thermal conductivity tensor scaling laws of TPLS-based cellular structures are added to the optimization algorithm. This allows the optimized relative density distributions of the cellular structures distributed in a component to possess an accurate physical meaning.

Overview

The proposed design methodology for FGCS optimization mainly consists of four steps, which is illustrated by an example of the torque transfer optimization problem in the flow chart in Fig. 1. In step 1, the modified 3D implicit functions of several triply periodic minimal surfaces are developed to parametrically design different types of cellular structures and generate spatially graded cellular structures (see Sec. 4). In step 2, the numerical homogenization method is used to calculate the elastic tensor and the thermal conductivity tensor of TPLS-based cellular structures with different densities. According to the results, the equation of density is established to demonstrate the equivalent elastic tensor and the thermal conductivity tensor of the specific cellular structures (see Sec. 5). In step 3, the penalization function in the original SIMP-based optimization method is replaced by the calculated homogenized model. The accurate volume fraction of each unit is obtained by figuring out the tensor scaling laws added to the topology optimization problem (see Sec. 6.1). In step 4, the relative density of the numerical results of structure optimization is mapped into the parametric TPLS-based implicit functions, which realize an optimum lightweight cellular structure (see Sec. 6.2). Finally, FEA simulation is conducted to verify the structural robustness of our optimized results in Sec. 7, before conclusions are drawn in Sec. 8.

TPLS-Based Cellular Structures Design

Mathematical Description.

Level surfaces are represented by function F : R3R at point (x, y, z) ∈ R3, which satisfies the equation
F(x,y,z)=t
(1)
where t is the level parametric. The form of the surface is controlled by F(x, y, z), while the parameter t determines the volume fraction of the region. A periodic surface F(x, y, z) can generally be defined as [43]
F(x,y,z)=hkl|ϕ(hkl)|cos(2πL(hx+ky+lz)αhkl)
(2)
where ϕ(hkl) denotes the structure magnitude factor reflecting the symmetry of the structure; the set of allowed hkl values express the reciprocal cellular vectors of a given lattice [44]; αhkl represents the phase shift, and L is the edge length of periods. Therefore, specific periodic structures and phases can be constructed based on the implicit form.
In this paper, four different cellular structures based on typical triply periodic minimal surfaces are investigated. The cellular types are the gyroid (G), diamond (D), primitive (P), and iWP (W), which have incremental nodal connectivity with threefold, fourfold, sixfold, and eightfold connected nodes, respectively, as shown in Fig. 2. TPLS represents the solution to scalar-valued functions of three independent variables. The three independent variables can be considered as the x, y, and z coordinates of a point in three-dimensional Euclidean space. The close approximation to the G, D, P, and W level surfaces [44] are composed of trigonometric functions F : R3R at point (x, y, z) ∈ R3, and they are defined as follows:
FG(x,y,z)=sin(x)cos(y)+sin(y)cos(z)+sin(z)cos(x)t
(3)
FD(x,y,z)=sin(x)sin(y)sin(z)+sin(x)cos(y)cos(z)+cos(x)sin(y)cos(z)+cos(x)cos(y)sin(z)t
(4)
FP(x,y,z)=cos(x)+cos(y)+cos(z)t
(5)
FW(x,y,z)=2(cos(x)cos(y)+cos(y)cos(z)+cos(z)cos(x))+(cos(2x)+cos(2y)+cos(2z))t
(6)
The level parameter t is a variable and determines the volume fractions pertaining to the regions separated by the surface. The cellular structures are subjected to translation and scaling. Moreover, they are pinched off into unconnected blobs as the value of the parameter t becomes sufficiently negative or positive. The relative density, or the volume fraction of the solid phase, is an important parameter for cellular materials. It is usually used as a variable to express other properties, such as elastic modulus and yield stress. It is essential to acquire the relationship between parameter t and relative density ρTPLS for TPLS-based cellular structures generated from Eqs. (3)(6), respectively. The data plotted in Fig. 3 reflect how the relative density of G, D, P, and W structures change with the parameter t. The best fitting of data gives the scaling laws of relative density for G, D, P, and W structures, and they are given as follows:
ρG=0.331×t+0.495,1.5t1.5
(7)
ρD=0.393×t+0.483,1.232t1.232
(8)
ρP=0.058×t2+0.269×t+0.575,1.513t3.12
(9)
ρW=0.029×t2+0.279×t+0.527,1.556t1.9
(10)

The most obvious trend in Fig. 3 is that the relative density increases gradually with an increase in parameter t. It indicates the linear relationships for both G- and D-based cellular structures and the quadratic polynomial relations for both P- and W-based cellular structures. According to the error analysis of the fitting, R2 values are very close to 1.0, so the fitting accuracy of Eqs. (7)(10) is very high.

Generating Graded Cellular Structures.

In this work, to generate TPLS-based continuous gradient cellular structures in 3D space, the level parameter t is modified as a function: t(x, y, z) in R3. Therefore, the TPLS-based density graded cellular can be expressed as follows:
FG(x,y,z)=sin(x)cos(y)+sin(y)cos(z)+sin(z)cos(x)t(x,y,z)
(11)
FD(x,y,z)=sin(x)sin(y)sin(z)+sin(x)cos(y)cos(z)+cos(x)sin(y)cos(z)+cos(x)cos(y)sin(z)t(x,y,z)
(12)
FP(x,y,z)=cos(x)+cos(y)+cos(z)t(x,y,z)
(13)
FW(x,y,z)=2(cos(x)cos(y)+cos(y)cos(z)+cos(z)cos(x))+(cos(2x)+cos(2y)+cos(2z))t(x,y,z)
(14)
In these functions, each value of t(x, y, z) determines the density of a certain point (x, y, z) ∈ R3, and it can be represented with an expression of high degree to achieve more complex density distribution. In this work, the calculated stress vector or temperature vector is transferred as a 3D matrix ρi,j,k to generate a complex density distribution. According to Eqs. (7)(10), the value of each parameter t corresponding to the density three-dimensional matrix ρi,j,k can be obtained, and then the parameter matrix ti,j,k can be established in 3D space. For example, Fig. 4 shows an example of a linear grading for G, D, P, and W unit cells where density varies from 0.1 at the right to 0.7 at the left of the cellular structures. In addition, material properties such as stress, Young’s modulus, and heat distribution are functions of their relative density in R3 [45]. Accordingly, if the mechanical property of a point in space is determined, the density value of the point can be obtained. The density values for all points can form a three-dimensional matrix, which corresponds to the values in t(x, y, z). Therefore, functional gradient structures that meet the mechanical or thermal properties can be generated by adding this 3D matrix into Eqs. (11)(14). In the next step, the developed parametric functionally graded TPLS-based cellular generation algorithm will be combined with structural optimization techniques.

Homogenized Model of TPLS-Based Cellular Structures

In this work, the effective properties of lattice material can be calculated by the asymptotic homogenization (AH) approach [46]. The AH method with a sound mathematical foundation has been widely used to calculate the effective properties of the periodic microlattice structure, the random porous structure, composite material, and also density-based topology optimization [46]. The AH approach is based on the analysis of a representative volume element of the cellular material at mesoscale. As shown in Fig. 5, it is the basic principle of the AH method. Suppose a body Ωε is infilled with a periodic cellular structure subjected to the traction t at the traction boundary Γt, a displacement d at the displacement boundary Γd, and a body force f is distributed in the domain. The cellular structure is considered to be a continuum solid domain Ω with equivalent properties when the AH method is used. Therefore, the computational cost is significantly reduced compared with the full-scale simulation of the explicitly represented cellular structure.

Effective Mechanical Elasticity Tensor.

In this work, the homogenized elasticity tensor are calculated by the homogenization theory [47] for TPLS-based unit cellular structures. The AH method assumes that each mesoscale unit cell in the macrostructure follows a periodic boundary conditions and outer forces corresponding to uniform strain fields is need to be analyzed. The measurable field quantity of a unit cellular structure u is the superposition of the macroscale quantity u0(x, y) and a small periodically fluctuated mesoscale quantity u1(x, y), which can be represented by the first-order asymptotic expansion [46]:
uε=u0(x,y)+εu1(x,y)+ε2u1(x,y)+y=x/ε,ε1
(15)
For linear elasticity, u0 is the average displacement on the macroscopic scale, and u1 and u2 represent the perturbation displacement on the microscale [46]. Homogenized elasticity tensor CTPLSH of a periodic unit cell can be rewritten in an equivalent discrete form in terms of elemental energy:
CTPLSH=1|V|e=1nVe(IBeχe)TDe(IBeχe)dVe
(16)
where the summation represents the assembly of n finite elements. The V is the unit cell volume, Ve is the volume of element e, I is a six times identity matrix, Be is the element strain-displacement matrix, χe is the matrix containing the element displacement vectors, and De is the constitutive matrix for the element, which for the basic solid material of the object. The global stiffness matrix K of the unit cellular microstructure can be described as follows:
K=e=1nVeBeTDeBedVe
(17)
The load fij that corresponds to macroscopic volumetric straining can be expressed as follows:
fij=e=1nVeBeTDeεijdVe
(18)
where the strains ɛij are chosen to be
ε11=(1,0,0,0,0,0)T,ε22=(0,1,0,0,0,0)T,ε33=(0,0,1,0,0,0)Tε12=(0,0,0,1,0,0)T,ε23=(0,0,0,0,1,0)T,ε13=(0,0,0,0,0,1)T
(19)
Finally, the displacement fields are computed by solving the finite element problem with six load cases by the following function:
[e=1nVeBeTDeBedVe]χij=e=1nVeBeTDeεijdV
(20)
where χij is the global displacement vector of the TPLS-based cellular structures. When the global displacements have been obtained, the homogenized constitutive matrix CijklH can be found as follows:
CijklH=1|V|e=1nVe(χe0(ij)χe(ij))TDe(χe0(kl)χe(kl))dVe
(21)
where χe0 contains the six displacement fields corresponding to the solid unit strains, and χe(ij) contains six columns corresponding to six displacement fields resulting from globally enforcing the unit strains in Eq. (20). The homogenized elasticity tensor is obtained based on periodic boundary conditions. For the TPLS-based cellular structures, this tensor is a symmetric matrix of the following form:
CTPLSH=[C1111HC1122HC1133HC1112HC1123HC1131HC2222HC2233HC2212HC2223HC2231HC3333HC3312HC3323HC3331HC1212HC1223HC1231HsymC2323HC2331HC3131H]
(22)

Effective Thermal Conductivity Tensor.

Similarly, according to the homogenization theory for elasticity tensor, the homogenized thermal conductivity tensor kTPLSH of a discretized periodic unit cell can be written as
kTPLSH=1|V|e=1nVe(IBetTe)Tke(IBetTe)dVe
(23)
where n is the number of finite elements of the discretized unit cellular structure, Ve is the volume of the finite element, I is a six times identity matrix, Bet is the element temperature gradient matrix, Te is the matrix containing the element nodal temperature vectors, and ke is the element thermal conductivity tensor. For the thermal finite element analysis method, the thermal stiffness matrix Kt of the unit structure can be expressed as follows:
Kt=e=1nVeBetTkeBetdVe
(24)
The nodal heat flux vector of the unit cell can be expressed as
fit=e=1nVeBetTketidVe
(25)
where temperature gradients ti are chosen to be
t1=(1,0,0)T,t2=(0,1,0)T,t3=(0,0,1)T
(26)
Finally, the global temperature vectors are computed by solving the finite element problem with three load cases:
[e=1nVeBetTkeBetdVe]Ti=e=1nVeBetTketidVe
(27)
where Ti is the global temperature vector of the TPLS-based cellular structures. When the temperature vectors have been obtained, the entries in the homogenized thermal constitutive kijH can be found as
kijH=1|V|e=1nVe(Te0(i)Te(i))Tket(Te0(j)Te(j))dVe
(28)
where Te0(i) contains the three temperature fields corresponding to the solid unit strains, and Te(i) contains three columns corresponding to three temperature fields resulting from globally temperature vectors in Eq. (24). Ve is the unit cell volume. ket is the thermal conductivity tensor for each element e. The tensor is a symmetric matrix of the form
kTPLSH=[k11Hk12Hk13Hk22Hk23Hsymk33H]
(29)

Scaling Laws of the Tensors.

Relative density ρr is an important characteristic of cellular structures. The “Gibson–Ashby” model [45] suggests that the relative elastic tensor of a cellular structure is a function of its relative density. This function is known as “scaling law.” This law can equally be applied to TPLS-based cellular structures [16,45]. Based on this theory, the relationship between the elastic or thermal tensors, and their relative density of the TPLS-based cellular structures is fitted, and this relationship is referred to as “tensor scaling laws” in this work.

In this paper, each coefficient CijklH and kijH can be expressed as a function of its relative density of the TPLS-based cellular structures. Because the TPLS-based cellular structures belong to cubic symmetric geometry, only C1111H, C1122H, C1212H, and k11H need to be calculated. To find the elastic and thermal conductivity tensor scaling laws of the TPLS-based cellular structures, different orders of polynomials and exponential functions are tried to fit the computational result to identify the relationship between constants of the elastic (thermal conductivity) tensors and the relative density of cellular structures. It is found that an exponential function not only guarantees that the function monotonically increases in the interval 0–1, but also the fitting error is small. Furthermore, it provides the best combination between accuracy and compactness of the scaling law. Therefore, the elastic tensor or the thermal conductivity tensor of cellular structures can be written as
CTPLSH(ρr)orkTPLSH=a1ea2ρr+a3
(30)
where ρr (ρr = ρc/ρ*, where ρc is the density of the cellular structure and ρ* is the density of the solid material) is the relative density of a cellular structure and a1, a2, and a3 are the constant symmetric matrices to be determined. For the anisotropic case, Eq. (30) gives the relationship of the mechanical properties and relative density. In the structure optimization procedure, the elastic tensor and the thermal conductivity tensor scaling laws become one when ρr is fully solid and becomes zero when ρr is void. Therefore, the fitted function passes two points: (1,1) and (0,0). Hence, the tensor scaling laws of cellular structures can be modified as
CTPLSH(ρr)orkTPLSH=a1ea2ρra1
(31)

The constant symmetric matrices a1 and a2 are obtained after fitting them to the homogenization calculation results. The scaling laws of the G, D, P, and W-based cellular structures are plotted in Fig. 6. The results show that the G-, D-, and W-based cellular structures have a relatively strong shear resistance, and the P-based cellular structure has a relatively strong pressure resistance. Therefore, the G, D, and W structures are suitable for optimized structural design with shear resistance requirements, while the P structure is suitable for optimized structural design with compressive requirements. In addition, the heat conduction performance of the G structure is better than that of other structures.

According to the calculated results, the elastic tensor and the thermal conductivity tensor scaling laws for the two constants of TPLS-based cellular structures are presented in Table 1, where C11111H, C1122H, C1212H, and k11H are the corresponding elastic constants and thermal conductivity constants of the composing material.

Structural Optimization

Objective Functions.

General topology optimization methods provide an optimal solution for retaining useful material within a prescribed design domain to achieve the required structural or thermal performances. The elements are considered as isotropic materials. However, almost all types of cellular structures belong to cubic symmetric geometry. The cubic symmetric geometric structure is a highly symmetrical geometrical structure. Furthermore, their mechanical properties are usually anisotropic and are symmetrical only in the direction of the three axes. Therefore, when the anisotropic cellular structure has a design material, Liu et al. [48] proposed a porous anisotropic material with penalization topology optimization method. Wang et al. [49] proposed a topological optimization with connectable graded microstructures, and the anisotropic properties of the microstructure are taken into account. Besides, Liu et al. [50], based on the MMC/MMV topology optimization framework, designed an orientation controllable lattice structure with an anisotropic mechanical performance. In this section, the minimum compliance cellular structure is obtained by adding the calculated elastic tensor or thermal conductivity tensor scaling laws into the optimization algorithm. Therefore, the structure optimization problem can be treated as the optimization of density distribution of the cellular structures of a mechanical part.

For the mechanical stiffness problem, an algorithm is proposed to optimize the infilling process by considering the materials properties of cellular structures in this work. The objective function used here is to identify the optimal distribution of material density, subjected to target volume constraint, to minimize structural compliance. Accordingly, the mathematic expression of the minimum compliance problem for the cellular structure utilizes the following form:
find:ρminc(ρ)=12UTKTPLS(ρ)Us.t.CTPLSH=CTPLSH(ρ)KTPLS(ρ)U=fi=1mviρeiV0<ρminρeiρmax1given:KTPLS(ρ)=e=1nVeBeTCTPLSH(ρ)BedVe
(32)
where structural compliance c(ρ) is the objective function, ρ is the vector of the element relative densities in the design domain and the element densities that constitute ρ are the design variables, U is the displacement vector, KTPLS(ρ) is the global stiffness matrix, and f is the external load vector. ρe is the relative density of element e, while vi is its corresponding element volume. The design variable ranges from ρmin to ρmax. The first constraint is the effective elastic tensor scaling law, which is established in Sec. 5. The second constraint is the Hooke’s law equilibrium equation. The third constraint limits the total design volume to V. The final constraint is the minimum and maximum design variables for the relative density. The optimization problem is solved by the optimality criteria algorithm in this work.
For the thermal conductivity problem, the goal is to find the optimal conduction path to remove heat from a design domain. This can also be seen as a problem where the minimization of the mean temperature in the design domain is sought. The complete optimization problem can be written by means of a finite element formulation as
find:ρminct(ρ)=12TTKTPLSt(ρ)Ts.t.kTPLSH=kTPLSH(ρ)KTPLSt(ρ)T=fti=1mviρeiV0<ρminρeiρmax1given:KTPLSt(ρ)=e=1nVeBetTkTPLSH(ρ)BetdVe
(33)
where ct(ρ) is the objective function, ρ is the vector of the element relative densities in the design domain, KTPLSt(ρ) is the global conductivity matrix, T is the global nodal temperature vector, and ft is the applied heat load vector. ρe is the relative density of element e, while vi is its corresponding element volume. The design variable ranges from ρmin to ρmax. The first constraint is the thermal conductivity tensor scaling law, which is established in Sec. 5.2. The second constraint is the Fourier’s law equilibrium equation. The third constraint is the limit design volume. The final constraint is the minimum and maximum design variables for the relative density.

Relative Density Mapping Method.

In this work, the RDM methods are used to match the mechanical properties of the optimized object. The RDM methods employ the corresponding density values of finite elements calculated from the structure optimization process. The mapping process is shown in Fig. 7. First, a design space XR3 is given with target volume constraint, boundary conditions, and the cellular structure xR3 used for infill. Second, the effective cellular structure mechanical properties P(X) can be calculated via map h. In this work, map h represents the numerical homogenization method. Third, the established scaling law of cellular structures is added in the optimization procedure, and this process is called mapping g. Therefore, the target density distribution of optimized material property fields P(Y) is achieved, and it is converted into a three-dimensional matrix. In the hmap1 mapping process, the density of each element in the optimized output is assigned to all nodes corresponding to that element. Then, the resulting density on each node is calculated by using Eq. (34) to get the average density for all elements of the shared pattern.
ρijk=XijkNnode
(34)
where ρijk is the synthetic density of the node point (i, j, k) and Xijk is the density of the elements that share the node point. Nnode is the number of elements sharing the node (Nnode can be 1, 2, 4, 8). Thus, Eq. (34) provides the average mapping density at each element node. After the density at the element node is obtained, the number of cells filled up in each direction can be controlled by linear interpolation across the entire element volume along three dimensions. For example, interpolation is performed between two points (xa, y, z) and (xb, y, z) adjacent in the x-axis direction, and the corresponding densities are ρxa and ρxb. Then, the distance between the two points is divided into n points {x1, x2, …, xn}, each point is (xi, y, z), and the density of each interpolation point is
ρxi=ρxb(ρxbρxa)xbxixbxai=1,2,,n
(35)

In the same way, arbitrary density interpolation can be performed throughout the design space to control the number of cells in a certain direction. Lastly, the matrix can be exported to the FGCS modeling method in Sec. 4.1, and the resulting structure Y can be generated. By combining the mapping h with the mapping g, a comprehensive mapping method f = hgh−1 is established. Therefore, the finally optimized results can be directly obtained.

Results and Discussion

In this work, the proposed gradient cellular structures design and optimization technique is programmed by the matlab code. To illustrate the effectiveness of the optimized gradient cellular structures herein, different design parameters are analyzed, comparison with existing typical gradient lattice design methods is made, and structural robustness is discussed, separately. The gyroid cellular structure is selected as a representative for the design and optimization in this section. For mechanical performance, it is assumed that the Young’s modulus for the solid material that constructs the cellular materials is 2.15 Gpa and the Poisson’s ratio is 0.3. For thermal conductive performance, for the solid material, the Young’s modulus of 60 GPa, the density of 2.5 × 10−6 kg/mm3, and the Poisson ratio of 0.35, the thermal conductivity k = 0.235 W/(mm K) are assumed.

Size-Effect Estimation.

The homogenization theory is based on a periodic distribution of cell structures that are much smaller than the macroscopic scale. The ideal condition is that the microscale unit size is much smaller than the macroscopic scale. However, when considering the results of actual optimization to meet the minimum size of existing AM, we need to find a minimum cell size that best satisfies the homogenization theory. In this paper, the influence of the size and boundary of the cellular structure on homogenization accuracy is discussed through numerical studies.

As shown in Fig. 8(a), a cubic block of 30 × 30 × 30 mm is employed as a benchmark with a nonperiodic boundary to study the influence of element size. The fixed boundary conditions are subjected to the bottom of the cube, while a load of 1000 N is applied to the upper surface while the other surfaces are unstressed, as shown in Fig. 8(b). For thermal conductivity, the temperature is 20° C at the top surface, a heat flux of 0.1 W/mm2 is applied to the bottom surface, and all other surfaces are perfectly insulated in Fig. 8(c). The cube is infilled with a gyroid-based cellular structure and the volume fraction is 0.3, and the constitutive model in Eq. (31) is used in the simulation. The analysis process is performed using the OptiStruct 2017 software. Accordingly, the maximum total deformation of the homogenized model is 0.092 mm and the maximum temperature of the homogenized model is 121.5° C.

To test the influence of the unit size on the nonperiodic boundaries in this benchmark, the number of unit cells in one direction is set as 1, 2, 3, 4, 5, 6, 8, 10, and 12 (Figs. 8(d)8(l)). Then, a full-scale simulation based on these models is performed to compare with the homogenized model. Figure 9(a) shows the maximum deformation of the upper surface relative to the number of unit cells in one direction. It can be seen from the results that the number of unit cells along one direction is larger than six; then the maximum total deformation gradually converges to the homogenized result. Besides, as can be seen in Fig. 9(b), as the unit cell size decreases, the temperature distribution of the bottom surface gradually approaches the homogenized result. In particular, when the number of unit cells in the heat transfer direction is larger than six, the difference between the homogenized model and the full-scale simulation is less than 2%. Therefore, six-unit cells along one direction of the design component are considered to be the minimum number in this work. In summary, the cell size and the relative density ranges have been appropriately determined to maximize the accuracy of the AH model while addressing manufacturability issues.

Structural Mechanical Performance

Example.

In this section, the classic 3D cantilever beam design problem is solved to demonstrate the effectiveness of the proposed methodology. Figure 10(a) shows that the dimensions of the beam model are 60 mm × 20 mm × 20 mm. The loading and constraints are also indicated in Fig. 10(a). The design domain is discretized using a 300 × 100 × 100 uniform grid. A mesh of 15,125 eight-node hexahedron elements is used to discretize the design domain and solve the FEA model. A load of 100 N is initially distributed over the object, as shown in the figure. The selected material parameters are the modulus of elasticity E = 2.15 GPa and Poisson’s ratio v = 0.3. The elastic tensor scaling law of the gyroid-based cellular structure applied to the optimization procedure is expressed in Eq. (30), and the correlation coefficients can be queried in Table 1. The target volume fractions are set to 0.3, 0.4, and 0.5, respectively. In addition, the minimum and maximum relative densities are set as uniform in Fig. 10(c), set as ρmin = 0.1 and ρmax = 0.9 in Fig. 10(e), and set as ρmin = 0.1 and ρmax = 0.8 in Fig. 10(f).

Figure 10(b) shows the optimized cellular structure by using gyroid to infill the design domain. The optimized structure has several notable characteristics. First, compared with the traditional topological optimization, the cellular structure is filled in the whole design area, which makes the overall structure robust, which will be discussed in the following sections. Second, the cellular structure is not only the endurance structure but also can be used as the support structure of the 3D printing process to improve manufacturability. Finally, it can be seen from the cross section area that the cellular structures with different densities can be smoothly connected (Fig. 10(c)), which means the smooth continuous structure can effectively transfer mechanical properties. In particular, the effects of different filling densities and different density intervals on the design results are also illustrated in Figs. 10(d)–10(f). The main conclusions are as follows: First, according to Figs. 10(e) and 10(f), the stiffness of the optimized structure with a density of 0.1–0.9 intervals is better than that of the density of 0.1–0.8 intervals. Second, by comparing Fig. 10(d) with Figs. 10(e) and 10(f), it can be seen that the optimized filling method makes the overall structure to have a better stiffness than that of uniform structures. For example, Fig. 11 illustrates the stiffness obtained by homogenization result and the full-scale simulation of Case #1 and Case #2. in Fig. 10. Table 2 lists the corresponding values of the stiffness of the linear elastic region for the nine cases in Fig. 10. From Fig. 11 and Table 2, it can be seen that the predicted stiffness value from the homogenization model is closer to the value of the full-scale simulation, and the stiffness of our optimized results is significantly better than that of the uniform form. This means that the proposed homogenization model can be effectively used to obtain a density variable lattice structure for optimizing the stiffness of an object.

Robustness.

The cellular structures possess an advantage of being able to tolerate localized damage. Parts in the manufacturing or working process are prone to local damage caused by stress concentration or an accidental event. These forms of damage mainly include accidental collision, fatigue failure, corrosion, and manufacturing defects. With traditional topology optimization, the structure retains the stressed material, and the rest is deleted. Therefore, when the optimized results are partially destroyed, the mechanical properties of the global structure will be seriously affected. The advantage of the density optimized cellular structures used in this work is that even if parts are locally broken, the overall stiffness can persist.

First, the half Messerschmitt-Bölkow-Blohm (MBB) beam (60 mm × 20 mm × 20 mm) is optimized by three different optimization strategies with the same materials and boundary conditions. They are the SIMP-based topology optimization method, the uniform gyroid-based cellular lightweight design method, and our proposed method, respectively (as shown in Fig. 12(a)). From the results, we can see that the stiffness of the topology optimization method is the best, and our method is slightly smaller than it. Moreover, compared with the uniform design method, the above two approaches have apparent advantages in stiffness. However, the SIMP-optimized results need much more support structures (vertical structures in Fig. 12(a)) for the AM process.

Second, to demonstrate the effectiveness of damage tolerance of distributed cellular structures, several simplified local deficiencies for the design results in the same place are added, as shown in Fig. 12(b). Subsequently, FEA simulations are completed separately for the original optimized models and the models with the inflicted flaws in different areas. The analysis process is performed using the OptiStruct 2017 software. The material properties used in the simulation are as follows: E = 2.15 GPa, v = 0.3, and F = 100 N. The maximum displacement of these models is counted in Fig. 12(c). For example, the displacement of the damaged structures (deficiency in A) and the nondamaged ones is compared; the displacement changes by 192.3% for the SIMP method, while it changes only by 33.6% for our optimized infill method. Furthermore, from the deformation results caused by different damaged areas, it can be seen that the displacement of the SIMP method fluctuates considerably. However, our proposed FGCS-based optimization method has less influence on its displacement under different failure modes. Although the displacement of a uniform structure under different failure modes is less affected, the overall deformation is always larger than that of other design strategies. The above result indicates that the proposed infill optimization technique produces a less damage-sensitive design than that of the standard SIMP method.

Comparison.

In this section, the proposed elastic tensor scaling law added optimization method will be compared with the cellular relative density mapping based on topology optimization with isotropic material [10,25]. In addition, it will be compared with the truss structure-based relative density optimization method [28,29]. The above methods correspond to “Method-1”, “Method-2”, and “Method-3” in Table 3, respectively. The first column of Table 3 is the design domain with the boundary conditions, as well as the items that need to be compared. The three models of Method-1, Method-2 and Method-3 are optimized with the same material properties, volume constraints and the compliance values for the two iterations are less than 0.01 as the same convergence condition. Besides, the maximum stress is obtained by FEA of the different optimized results.

According to the compliance values, it indicates that our proposed method is 1.35 times stiffer than Method-2 and 1.12 times stiffer than Method-3. Moreover, from the results of the finite element analysis, it can be seen that our proposed method has a maximum von Mises stress that is 52.6% of the maximum stress in Method-2 and 32.9% of the maximum stress in Method-3. This occurs because there is no direct material property relationship between the optimization calculation process in Method-2 and the filled cellular structure. In addition, the geometrical shape in Method-1 is smoother and more continuous than Method-3. Therefore, Method-1 is not an easy-to-produce stress concentration phenomenon. The above results show that the proposed method has better structural performance.

Thermal Conductive Performance.

Heat sinks are passive heat exchangers that are used widely in electronic devices to cool key components such as central processing units or graphic cards in computers. In this section, several heat-sink designs are presented. Their thermal conductive behaviors are compared. The initial geometry and boundary conditions are given in Fig. 13(a). It consists of a design domain 180 × 180 × 30 mm and has a sink on the middle of the top surface. A uniform surface heat flux of magnitude 0.1 W/mm2 is applied on the top surface and a forced convection exists over the structural boundaries. The optimization goal of the heat-sink is to minimize the thermal compliance satisfy a volume constraint of 30% of the design domain. For the solid material, the Young’s modulus of 60 GPa, the density of 2.5 × 10−6 kg/mm3, and the Poisson ratio of 0.35 are assumed. Other parameters are set as follows: thermal conductivity k = 0.235 W/(mm K) and prescribed temperature T0 = 20° C on the bottom, left, and right boundaries. Figure 13(b) illustrates a reference heat-sink design with uniform cellular structures and the corresponding thermal conductive simulation result. It has a maximum temperature of 66.1 °C. Figs. 13(c) and 13(d) show the cellular relative density distribution optimized results by our proposed method. The minimum and maximum relative densities are set as ρmin = 0.1 and ρmax = 0.8 in Fig. 13(c) and as ρmin = 0.1 and ρmax = 0.9 in Fig. 13(d). It can be seen from the simulation results that the optimized results are superior to the thermal conductivity of the homogeneous structure. Besides, in the optimal design results of the same volume fraction, the density range of 0.1 to 0.9 has better heat conduction performance than the density range of 0.1–0.8. Furthermore, the results of the homogenization method are very close to the full-scale simulation results in Table 4. It implies that the proposed homogenization model can be effectively used to obtain a density variable lattice structure for optimizing the temperature distribution.

Conclusions

In this work, a generative design algorithm combining the TPLS-based graded cellular structure modeling method with structure optimization is developed to optimize the parts for additive manufacturing. In contrast to the traditional uniform cellular design and standard topology optimization method, our technique not only optimizes the stiffness of the structure, but also makes the optimized result structural robustness concerning damage tolerant. Besides, the thermal conductivity of the cellular structures can also be effectively optimized by our method. The conclusions of this paper are as follows:

  • We presented an FGCS modeling methodology based on the TPLS-based cellular structures. It can be used to map material properties into 3D space.

  • A structure optimization strategy was established to achieve the optimum density distribution of cellular structures. The scaling laws of the elastic tensor (thermal conductivity tensor) of the TPLS-based cellular structures were calculated by using the numerical homogenization method, and it was added to the optimization process. This allows the optimized relative density distributions of the cellular structures distributed in a component to possess an accurate physical meaning.

In future work, we will expand our approach to the study of multiple physical coupling optimizations by using gradient cellular structures.

Acknowledgment

This work is supported by financial support from the Program of China Scholarships Council (No. 201706830039); National Science Foundation of China (Grant No. 51775273); Natural Science Foundation of Jiangsu Province, China (No. BK20161487); and Six talent peaks project in Jiangsu Province, China (Grant No. GDZB-034).

References

1.
Thompson
,
M. K.
,
Moroni
,
G.
,
Vaneker
,
T.
,
Fadel
,
G.
,
Campbell
,
R. I.
,
Gibson
,
I.
,
Bernard
,
A.
,
Schulz
,
J.
,
Graf
,
P.
,
Ahuja
,
B.
, and
Martina
,
F.
,
2016
, “
Design for Additive Manufacturing: Trends, Opportunities, Considerations, and Constraints
,”
CIRP Ann.-Manuf. Technol.
,
65
(
2
), pp.
737
760
.
2.
Orme
,
M. E.
,
Gschweitl
,
M.
,
Ferrari
,
M.
,
Madera
,
I.
, and
Mouriaux
,
F.
,
2017
, “
Designing for Additive Manufacturing: Lightweighting Through Topology Optimization Enables Lunar Spacecraft
,”
J. Mech. Des.
,
139
(
10
),
100905
.
3.
Dong
,
G.
,
Tang
,
Y.
, and
Zhao
,
Y. F.
,
2017
, “
A Survey of Modeling of Lattice Structures Fabricated by Additive Manufacturing
,”
J. Mech. Des.
,
139
(
10
),
100906
.
4.
Cheng
,
L.
,
Liu
,
J.
,
Liang
,
X.
, and
To
,
A. C.
,
2018
, “
Coupling Lattice Structure Topology Optimization With Design-Dependent Feature Evolution for Additive Manufactured Heat Conduction Design
,”
Comput. Methods Appl. Mech. Eng.
,
332
, pp.
408
439
.
5.
Cheng
,
L.
,
Liu
,
J.
, and
To
,
A. C.
,
2018
, “
Concurrent Lattice Infill With Feature Evolution Optimization for Additive Manufactured Heat Conduction Design
,”
Struct. Multidiscip. Optim.
,
58
, pp.
511
535
.
6.
Bates
,
S. R. G.
,
Farrow
,
I. R.
, and
Trask
,
R. S.
,
2016
, “
3D Printed Elastic Honeycombs With Graded Density for Tailorable Energy Absorption
,”
Active and Passive Smart Structures and Integrated Systems 2016
,
International Society for Optics and Photonics
, Vol.
9799
,
979907
.
7.
Schaedler
,
T. A.
, and
Carter
,
W. B.
,
2016
, “
Architected Cellular Materials
,”
Annu. Rev. Mater. Res.
,
46
, pp.
187
210
.
8.
Lu
,
L.
,
Sharf
,
A.
,
Zhao
,
H.
,
Wei
,
Y.
,
Fan
,
Q.
,
Chen
,
X.
,
Savoye
,
Y.
,
Tu
,
C.
,
Cohen-Or
,
D.
, and
Chen
,
B.
,
2014
, “
Build-to-Last: Strength to Weight 3D Printed Objects
,”
ACM Trans. Graph. (TOG)
,
33
(
4
),
97
.
9.
Li
,
D.
,
Dai
,
N.
,
Jiang
,
X.
, and
Chen
,
X.
,
2016
, “
Interior Structural Optimization Based on the Density-Variable Shape Modeling of 3D Printed Objects
,”
Int. J. Adv. Manuf. Technol.
,
83
(
9–12
), pp.
1627
1635
.
10.
Brackett
,
D.
,
Ashcroft
,
I.
, and
Hague
,
R.
,
2011
, “
Topology Optimization for Additive Manufacturing
,”
Proceedings of the Solid Freeform Fabrication Symposium
, Vol.
1
, pp.
348
362
.
11.
Aremu
,
A. O.
,
Maskery
,
I. A.
,
Tuck
,
C. J.
,
Ashcroft
,
I. A.
,
Wildman
,
R. D.
, and
Hague
,
R. J. M.
,
2016
, “
Effects of Net and Solid Skins on Self-Supporting Lattice Structures
,”
Challenges in Mechanics of Time Dependent Materials
, Vol.
2
, pp.
83
89
.
12.
Steven
,
G. P.
,
1997
, “
Homogenization of Multicomponent Composite Orthotropic Materials Using FEA
,”
Int. J. Numer. Methods Biomed. Eng.
,
13
(
7
), pp.
517
531
.
13.
Xie
,
Y. M.
,
Zuo
,
Z. H.
,
Huang
,
X.
, and
Yang
,
X.
,
2013
, “
Predicting the Effective Stiffness of Cellular and Composite Materials With Self-Similar Hierarchical Microstructures
,”
J. Mech. Mater. Struct.
,
8
(
5
), pp.
341
357
.
14.
Lord
,
E. A.
, and
Mackay
,
A. L.
,
2003
, “
Periodic Minimal Surfaces of Cubic Symmetry
,”
Curr. Sci.
,
85
, pp.
346
362
.
15.
Wu
,
X.
,
Ma
,
D.
,
Eisenlohr
,
P.
,
Raabe
,
D.
, and
Fabritius
,
H.-O.
,
2016
, “
From Insect Scales to Sensor Design: Modelling the Mechanochromic Properties of Bicontinuous Cubic Structures
,”
Bioinspir. Biomim.
,
11
(
4
),
045001
.
16.
Khaderi
,
S. N.
,
Deshpande
,
V. S.
, and
Fleck
,
N. A.
,
2014
, “
The Stiffness and Strength of the Gyroid Lattice
,”
Int. J. Solids Struct.
,
51
(
23–24
), pp.
3866
3877
.
17.
Michielsen
,
K.
, and
Stavenga
,
D. G.
,
2008
, “
Gyroid Cuticular Structures in Butterfly Wing Scales: Biological Photonic Crystals
,”
J. R. Soc. Interface
,
5
(
18
), pp.
85
94
.
18.
Gandy
,
P. J. F.
,
Cvijović
,
D.
,
Mackay
,
A. L.
, and
Klinowski
,
J.
,
1999
, “
Exact Computation of the Triply Periodic D (‘Diamond’) Minimal Surface
,”
Chem. Phys. Lett.
,
314
(
5–6
), pp.
543
551
.
19.
Rajagopalan
,
S.
, and
Robb
,
R. A.
,
2006
, “
Schwarz Meets Schwann: Design and Fabrication of Biomorphic and Durataxic Tissue Engineering Scaffolds
,”
Med. Image Anal.
,
10
(
5
), pp.
693
712
.
20.
Jung
,
Y.
, and
Torquato
,
S.
,
2005
, “
Fluid Permeabilities of Triply Periodic Minimal Surfaces
,”
Phys. Rev. E
,
72
(
5
),
056319
.
21.
Abueidda
,
D. W.
,
Al-Rub
,
R. K. A.
,
Dalaq
,
A. S.
,
Lee
,
D.-W.
,
Khan
,
K. A.
, and
Jasiuk
,
I.
,
2016
, “
Effective Conductivities and Elastic Moduli of Novel Foams With Triply Periodic Minimal Surfaces
,”
Mech. Mater.
,
95
, pp.
102
115
.
22.
Mahmoud
,
D.
, and
Elbestawi
,
M. A.
,
2017
, “
Lattice Structures and Functionally Graded Materials Applications in Additive Manufacturing of Orthopedic Implants: A Review
,”
J. Manuf. Mater. Process.
,
1
(
2
), p.
13
.
23.
Torres-Sanchez
,
C.
, and
Corney
,
J. R.
,
2009
, “
Toward Functionally Graded Cellular Microstructures
,”
J. Mech. Des.
,
131
(
9
),
091011
.
24.
Han
,
Y.
, and
Lu
,
W. F.
,
2018
, “
A Novel Design Method for Nonuniform Lattice Structures Based on Topology Optimization
,”
J. Mech. Des.
,
140
,
091403
.
25.
Panesar
,
A.
,
Abdi
,
M.
,
Hickman
,
D.
, and
Ashcroft
,
I.
,
2018
, “
Strategies for Functionally Graded Lattice Structures Derived using Topology Optimisation for Additive Manufacturing
,”
Addit. Manuf.
,
19
, pp.
81
94
.
26.
Liu
,
X.
, and
Shapiro
,
V.
,
2016
, “
Sample-Based Design of Functionally Graded Material Structures
,”
ASME 2016 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference
,
V02AT03A035
.
27.
Wu
,
J.
,
Wang
,
C. C. L.
,
Zhang
,
X.
,
Zhang
,
X.
, and
Westermann
,
R.
,
2016
, “
Self-Supporting Rhombic Infill Structures for Additive Manufacturing
,”
Comput. Aided Des.
,
80
, pp.
32
42
.
28.
Tang
,
Y.
,
Kurtz
,
A.
, and
Zhao
,
Y. F.
,
2015
, “
Bidirectional Evolutionary Structural Optimization (BESO) Based Design Method for Lattice Structure to Be Fabricated by Additive Manufacturing
,”
Comput. Aided Des.
,
69
, pp.
91
101
.
29.
Alzahrani
,
M.
,
Choi
,
S. K.
, and
Rosen
,
D. W.
,
2015
, “
Design of Truss-Like Cellular Structures Using Relative Density Mapping Method
,”
Mater. Des.
,
85
, pp.
349
360
.
30.
Daynes
,
S.
,
Feih
,
S.
,
Lu
,
W. F.
, and
Wei
,
J.
,
2017
, “
Optimisation of Functionally Graded Lattice Structures Using Isostatic Lines
,”
Mater. Des.
,
127
,
215
223
.
31.
Cheng
,
L.
,
Zhang
,
P.
,
Bai
,
J.
,
Robbins
,
J.
, and
To
,
A.
,
2017
, “
Efficient Design Optimization of Variable-Density Cellular Structures for Additive Manufacturing: Theory and Experimental Validation
,”
Rapid Prototyping J.
,
23
(
4
), pp.
660
677
.
32.
Wang
,
X.
,
Zhang
,
P.
,
Ludwick
,
S.
,
Belski
,
E.
, and
To
,
A. C.
,
2018
, “
Natural Frequency Optimization of 3D Printed Variable-Density Honeycomb Structure via a Homogenization-Based Approach
,”
Addit. Manuf.
,
20
,
189
198
.
33.
Bendsøe
,
M. P.
, and
Kikuchi
,
N.
,
1988
, “
Generating Optimal Topologies in Structural Design Using a Homogenization Method
,”
Comput. Methods Appl. Mech. Eng.
,
71
(
2
), pp.
197
224
.
34.
Bendsøe
,
M. P.
,
1989
, “
Optimal Shape Design as a Material Distribution Problem
,”
Struct. Optim.
,
1
(
4
), pp.
193
202
.
35.
Huang
,
X.
, and
Xie
,
Y. M.
,
2007
, “
Convergent and Mesh-Independent Solutions for the Bi-Directional Evolutionary Structural Optimization Method
,”
Finite Elem. Anal. Des.
,
43
(
14
), pp.
1039
1049
.
36.
Wang
,
M. Y.
,
Wang
,
X.
, and
Guo
,
D.
,
2003
, “
A Level Set Method for Structural Topology Optimization
,”
Comput. Methods Appl. Mech. Eng.
,
192
(
1–2
), pp.
227
246
.
37.
Guo
,
X.
,
Zhang
,
W.
, and
Zhong
,
W.
,
2014
, “
Doing Topology Optimization Explicitly and Geometrically—A New Moving Morphable Components Based Framework
,”
J. Appl. Mech.
,
81
(
8
),
081009
.
38.
Guo
,
X.
,
Zhang
,
W.
,
Zhang
,
J.
, and
Yuan
,
J.
2016
, “
Explicit Structural Topology Optimization Based on Moving Morphable Components (MMC) With Curved Skeletons
,”
Comput. Methods Appl. Mech. Eng.
,
310
, pp.
711
748
.
39.
Zhang
,
W.
,
Yang
,
W.
,
Zhou
,
J.
,
Li
,
D.
, and
Guo
,
X.
,
2017
, “
Structural Topology Optimization Through Explicit Boundary Evolution
,”
J. Appl. Mech.
,
84
(
1
),
011011
.
40.
Zhu
,
Y.
,
Li
,
S.
,
Du
,
Z.
,
Liu
,
C.
,
Guo
,
X.
, and
Zhang
,
W.
,
2019
, “
A Novel Asymptotic-Analysis-Based Homogenisation Approach Towards Fast Design of Infill Graded Microstructures
,”
J. Mech. Phys. Solids
,
124
, pp.
612
633
.
41.
Xia
,
L.
,
Xia
,
Q.
,
Huang
,
X.
, and
Xie
,
Y. M.
,
2016
, “
Bi-Directional Evolutionary Structural Optimization on Advanced Structures and Materials: A Comprehensive Review
,”
Arch. Comput. Methods Eng.
,
25
, pp.
437
478
.
42.
Sigmund
,
O.
, and
Maute
,
K.
,
2013
, “
Topology Optimization Approaches
,”
Struct. Multidiscip. Optim.
,
48
(
6
), pp.
1031
1055
.
43.
Von Schnering
,
H. G.
, and
Nesper
,
R.
,
1991
, “
Nodal Surfaces of Fourier Series: Fundamental Invariants of Structured Matter
,”
Z. Phys. B Condens. Matter
,
83
(
3
), pp.
407
412
.
44.
Gandy
,
P. J. F.
,
Bardhan
,
S.
,
Mackay
,
A. L.
, and
Klinowski
,
J.
,
2001
, “
Nodal Surface Approximations to the P, G, D and I-WP Triply Periodic Minimal Surfaces
,”
Chem. Phys. Lett.
,
336
(
3–4
), pp.
187
195
.
45.
Ashby
,
M. F.
,
2006
, “
The Properties of Foams and Lattices
,”
Philos. Trans. R. Soc. Lond. A Math. Phys. Eng. Sci.
,
364
(
1838
), pp.
15
30
.
46.
Bensoussan
,
A.
,
Lions
,
J. L.
, and
Papanicolaou
,
G.
,
2011
,
Asymptotic Analysis for Periodic Structures
,
American Mathematical Society
,
Providence, Rhode Island
.
47.
Andreassen
,
E.
, and
Andreasen
,
C. S.
,
2014
, “
How to Determine Composite Material Properties Using Numerical Homogenization
,”
Comput. Mater. Sci.
,
83
, pp.
488
495
.
48.
Liu
,
L.
,
Yan
,
J.
, and
Cheng
,
G.
,
2008
, “
Optimum Structure With Homogeneous Optimum Truss-Like Material
,”
Comput. Struct.
,
86
(
13–14
), pp.
1417
1425
.
49.
Wang
,
Y.
,
Chen
,
F.
, and
Wang
,
M. Y.
,
2017
, “
Concurrent Design With Connectable Graded Microstructures
,”
Comput. Methods Appl. Mech. Eng.
,
317
, pp.
84
101
.
50.
Liu
,
C.
,
Du
,
Z.
,
Zhang
,
W.
,
Zhu
,
Y.
, and
Guo
,
X.
,
2017
, “
Additive Manufacturing-Oriented Design of Graded Lattice Structures Through Explicit Topology Optimization
,”
J. Appl. Mech.
,
84
(
8
),
081008
.