Self-supporting Topology Optimization for Additive Manufacturing
Dengyang Zhao Ming Li Yusheng Liu












 































𝐱
𝐛̱
𝐟
Ç
𝒮
𝒰
ℳ
ℋ̋

e
b
1/2
𝐭
𝐧
𝐮̆
𝐯̌
𝐈
𝐕

𝐊
ℝ

𝐀
𝐆
𝐅
S
ℤ
P

τ̅



    thmTheoremlemLemmaelsarticle-num
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    cadaddress]Dengyang Zhaocadaddress]Ming Limycorrespondingauthorcadaddress]Yusheng Liu[mycorrespondingauthor]Corresponding author: liming@cad.zju.edu.cn[cadaddress]State Key Laboratory of CAD&CG, Zhejiang University, Hangzhou, China
    The paper presents a topology optimization approach that designs an optimal structure, called a self-supporting structure, which is ready to be fabricated via additive manufacturing without the usage of additional support structures. Such supports in general have to be created during the fabricating process so that the primary object can be manufactured layer by layer without collapse, which is very time-consuming and waste of material.
    
    The proposed approach resolves this problem by formulating the self-supporting requirements as a novel explicit quadratic continuous constraint in the topology optimization problem, or specifically, requiring the number of unsupported elements (in terms of the sum of squares of their densities) to be zero. Benefiting form such novel formulations, computing sensitivity of the self-supporting constraint with respect to the design density is straightforward, which otherwise would require lots of research efforts in general topology optimization studies. The derived sensitivity for each element is only linearly dependent on its sole density, which, different from previous layer-based sensitivities, consequently allows for a parallel implementation and possible higher convergence rate. In addition, a discrete convolution operator is also designed to detect the unsupported elements as involved in each step of optimization iteration, and improves the detection process 100 times as compared with simply enumerating these elements. The approach works for cases of general overhang angle, or general domain, and produces an optimized structures, and their associated optimal compliance, very close to that of the reference structure obtained without considering the self-supporting constraint, as demonstrated by extensive 2D and 3D benchmark examples.
    
    self-supporting  topology optimization explicit quadratic constraints  additive manufacturing  discrete convolution
    
    
    
    
    
    
    § INTRODUCTION
    
    Topology optimization aims to generate an optimal material distribution within a design domain under certain geometric or physical constraints. Since its introduction in late the 1980s <cit.>, this problem has attracted wide industrial and academic interests due to its large potentiality in engineering applications and its intrinsic mathematical challenges. Topology optimization has developed in many different forms, such as: homogenization <cit.>, density (SIMP) <cit.>, evolutionary approaches (BESO) <cit.>, level set <cit.>, or more recently IGA (iso-geometric analysis) <cit.>, to name a few. See <cit.> for a recent and comprehensive review on this topic.
    
    The complex geometric designs produced by topology optimization show the approach's superiority in balancing the geometric distribution and the target physical performance. Such designs are however very difficult to be manufactured directly via traditional subtractive or formative manufacturing techniques <cit.>. On the other hand, rapidly developing additive manufacturing technologies have the promise to overcome the barrier between the potentiality that the topology optimization approaches can provide and the limitations that traditional manufacturing technologies can fabricate. In reality, additive manufacturing is a natural counterpart to topology optimization in that they have very versatile capability to quickly generate and realize new components not existing before <cit.>. 
    
    Despite the enhanced geometric freedom associated with additive manufacturing, specific design rules must still be satisfied in order to ensure manufacturability. The fabrication overhang angle is an example of a rule which is of paramount importance to ensure that the part will not collapse when fabricating the designed structure layer by layer. A structure satisfying such an overhang angle constraint is called self-supporting. For example, Thomas <cit.> identified 45 degree as the typical maximum overhang angle with a large number of experiments. For a non self-supporting structure, its geometry has to be modified or additional support structures need to be generated. Modifying the geometry will ultimately reduce the structure's physical performance, while additional support raises the issue of automatic and minimum volume support design <cit.>, and further post-processing activities to remove the unwanted supports. In the case that the support is made of the same material as the main component, such as the selective laser melting (SLM) process using metals, it is extremely difficult to remove out the support structure. Particularly, when the generated supports are embedded within a closed volume of the model, it is impossible to remove them.
    
    The best strategy to resolve the issue of topology optimization for additive manufacturing is perhaps to design a completely self-supporting structure, via topology optimization, that can be fabricated directly without the usage of support materials. Brackett et al first suggested including the overhang angle constraints into the topology optimization process <cit.>, but does not produce a complete self-supporting structure. The first self-supporting structure built from topology optimization is due to the pioneering work of Gaynor and Guest in 2014 <cit.> (and a very recent journal version <cit.>), which is achieved via introducing a wedge-shaped filter during the topology optimization process. Also very recently, an excellent work was also conducted by Langelaar in 2D <cit.> and in 3D <cit.> via introducing a novel self-supporting filter into the topology optimization process, which is achieved via building smooth approximation to the minimum and maximum functions. Impressive 2D and 3D examples were also shown in these studies <cit.>.
    
    In this paper, an alternative novel self-supporting topology optimization approach is proposed to generate a structure of optimal physical performance that does not need any additional support materials. It is achieved via carefully formulating the self-supporting constraint as an explicit quadratic function with respect to the design density, specifically, requiring the number of unsupported elements (in terms of the sum of their densities) to be zero. Benefiting from the novel quadratic formulation, the self-supporting sensitivity for each element is straightforward to compute, and is only linearly dependent on density of the element itself; notice that designing a proper filter and computing the associated sensitivity usually requires lots of research efforts in general topology optimization framework <cit.>. The derived sensitivity does not involve density information of any other elements, and thus allows for a parallel implementation, which is particularly important for 3D problem of high DOFs. Previous approaches <cit.> have a nonlinear layer-based sensitivity expressions, and may thus inhibit parallel implementations, as also explained by the authors <cit.>. In addition, a discrete convolution operator is also designed to detect the unsupported elements as involved in each step of optimization iteration, and improves the detection process 100 times as compared with simply enumerating these elements.
    
    
    Comparisons between the proposed approach and previous studies <cit.> are also summarized in Table 1.
    
    
    
    The remainder of the paper is organized as follows. The novel formulation of self-supporting topology optimization, together with its overall numerical procedure, is presented in Section <ref>. Several numerical techniques behind the proposed approach are detailed in Section <ref>. Extensive 2D and 3D examples are demonstrated in Section <ref>. pFinally the paper is concluded in Section <ref>.
    
    
    
    § PROBLEM STATEMENT AND APPROACH OVERVIEW
    
    In the section, the self-supporting constraint is formulated as a quadratic continuous function in terms of the element density, and is integrated within the classical SIMP framework <cit.> for self-supporting topology optimization. Following on from this, the proposed numerical approach to resolve the problem is outlined.
    
    
    
     §.§ Supported and unsupported elements
    
    The supported elements generally stand for the structural elements that can be fabricated via an additive manufacturing technology without collapse with respect to the fabrication process. They are defined here using the concepts of a maximum printable supporting angle, or overhang angle, which is first assumed to be 45 degree following previous study <cit.>. Extensions of the approach to general overhang angles are also explained later. We also assume that the printing direction is following the positive-y axis direction in both 2D and 3D for ease of explanation.
    
    
    
    First consider a 2D discrete structured mesh model  consisting of square elements e(n,m), that is,
    
        ={e(n,m)| 1≤ n≤ N, 1≤ m≤ M },
    
    where n,m are the indices increasing along the x and y axes respectively. Without confusion, we also use e to represent a square element without explicitly mentioning its indices n,m. In addition, a density matrix  of size N× M is also associated to , where an entry value ρ(n,m)=1 or 0 respectively represents a solid or void element e(n,m) of .
    
    As illustrated in Fig. <ref>(a), given a solid element e(n,m) in  (in orange), it is supported, or called a supported element, if one of the three blue elements below it is solid. We formulate the self-supporting condition in a continuous form as follows: an element e(n,m)∈ is supported if
    
        ∑_n-1≤ r≤ n+1ρ(r,m-1)> 0.
    
    
    Correspondingly, the supporting set_S of model  is the set of all supported elements within , that is,
    
        _S={e(n,m)∈ | m=1  ∑_n-1≤ r≤ n+1ρ(r,m-1)> 0}.
    
    
    Similarly, given a 3D structured mesh model  consisting of cubic elements,
    
        ={e(n,m,l)| 1≤ n≤ N, 1≤ m≤ M, 1≤ l≤ L},
    
    where n,m,l are the indices increasing along the x,y,z axes respectively, the supporting set of  is similarly defined (see also Fig. <ref>(b)):
    
        _S= {e(n,m,l)∈ | m=1 ∑_n-1≤ r≤ n+1,
         l-1≤ s≤ l+1ρ(r,m-1,s)> 0}.
     
    Note also here that only five elements are included here as other element do not form an appropriate overhang angle with the orange element.
    
    Correspondingly, the set of unsupported elements of a model  is
    
        _U=∖_S.
    
    
     §.§ Formulation of self-supporting topology optimization
    
    The self-supporting topology optimization problem aims to find the optimal material distribution within a design domain under certain boundary conditions. As widely studied before, the problem of minimum compliance or equivalently maximum stiffness is examined here. Following the classical SIMP framework <cit.>, the problem of self-supporting topology optimization is formulated here as an optimization problem with an additional explicit self-supporting constraint. The constraint is carefully reformulated using a simple quadratic function with respect to the density, specifically, requiring the number of unsupported elements (in terms of the sum of square of their densities) is zero. Details are explained below.
    
    The problem of self-supporting topology optimization is stated as: find the optimal density distribution ,
    
        min_∈^N× Mc (,̆),  s.t.
    
        {[                                     ( )=̆( ),;                                  V()/V_0 ≤ f,;                        U() = ∑_e ∈_Uρ_e^2 ≤ϵ,;                         0< ρ_e≤ 1,  e=1,…, s, ].
    
    where  is the vector of design variables (element densities) to be computed, $̆ is the vector of global displacements andis the global stiffness matrix. The objective functionc(,̆)is the structure's compliance, defined as
    
        c(,̆)=^̆T.̆()is the nodal force vector,V()andV_0are the material volume and design domain volume,fis the prescribed volume fraction,_Uis the index set of unsupported elements as defined in (<ref>), andϵ>0is a small parameter closed to0. A penalty parameterp, usually set asp=3, is applied here for the 0,1 convergence of, or specifically,
    
        _e=ρ_e^p_e^0,
    
    where_e^0is the element stiffness matrix associated with an elementein the modelandρ_ethe associated element density.
    
    The only difference of the above conventions in (<ref>) with previous SIMP-based formulations is that it has an additional constraintU() = ∑_e ∈ρ_e^2 ≤ϵto meet the self-supporting requirement. This condition is based on the observation that when the sum of the element densities of unsupported elements tends to0, all the elements are self-supported; the square is used here so that its derivative is not constant. The simple quadratic expression allows for a straightforward sensitivity derivation of the self-supporting constraint, and ultimately results in a linear expression. 
    
    
    
     §.§ Approach overview
    
    
    The self-supporting topology optimization problem (<ref>) is ready-to-solve using the MMA approach noticing that the sensitivity is straightforward to compute (as can be further seen in Section <ref>). On the other hand, in order to further improve the approach's convergence and computational efficiency, the overall optimization process is carefully designed, as plotted in Fig. <ref> and detailed below.
    
    Firstly, whether an element is self-supporting is dependent on the printing direction. Different printing directions produce different optimization structures. Thus, if the printing direction is chosen arbitrarily, it may produce a structure totally different from the benchmark support-needed structure obtained without considering the self-supporting constraint, with a possible worse physical performance. In some very special case, the optimization approach may not converge. Thus, an appropriate printing direction is first set via generating a coarse structure via topology optimization without considering the self-supporting constraints. It chooses the coordinate axis direction with the least number of unsupported elements as the print direction.
    
    The criteria to generate this initial coarse structure, or to stop the initial optimization process stops, is based on the measure of non-discretenessM_nd, proposed in <cit.>. It stops when
    
        M_nd<M_nd^0,
    
    where
    
        M_nd=∑_e=1^N4ρ_e(1-ρ_e)/n× 100%,
    
    andnis the number of elements of the domain,ρ_eis the density an elemente. The valve parameterM_nd^0is set0.36, which corresponds to a structure of average density of0.1or0.9.
    
    The generated coarse structure is then set as an initial structure for self-supporting topology optimization, together with the chosen printing direction. During the optimization process, in order to prevent the self-supporting constraint hindering the formation of load-carrying structures, the constraint is imposed in a soft way using a similar strategy as previously performed in work <cit.>. Specifically, the value of the tolerance parameterϵin (<ref>) is decreased continuously from a relatively large value to the predefined toleranceϵ_0during the course of optimization.
    
    During the process of self-supporting topology optimization, a black-white filter is additionally applied in the last few steps to produce a totally 0-1 density distribution so that the discrete convolution finds exactly the self-supporting elements without the influence of gray elements. This is achieved using the Heaviside projection filter as originally designed by Guest et al <cit.>. Otherwise, gray elements may not be strong enough to support a black element above it, and special care has to be considered <cit.>.
    
    Lastly, in order to have a completely self-supporting structure (without any unsupported elements), the self-supporting constraint is added in a strict way in the last few iteration steps (when the number of unsupported elements is no longer decreasing) via removing non-self supporting elements (whose number is at most 5 is all tests given in this paper). This strategy has an ignorable influence on the final structure's compliance, noticing the nature of optimization approach, i.e. the approach may fluctuate between self-supporting constraints and target optimization and the volume faction constraint. Consequently, all these remaining elements are not essential in determining the structure's physical compliance. This is very different from removing unsupported elements from the support-needed structure in a post-processing step.
    
    
    
    
    § NUMERICAL ASPECTS
    
    
     §.§ Sensitivity analysis
    
    The sensitivity of the self-support constraintU()is straightforward to compute from (<ref>), and given as:
    
        ∂ U/∂ρ_e=
            
            2ρ_e   if  e ∈_U 
        
            0   if  e ∉_U,
    
    where_Uis the set of unsupported elements. 
    
    Note here that the sensitivity for each elementeis only dependent on the densityρ_eof the element itself, without information of density of any other elements, and can be easily implemented in parallel. Previous approaches <cit.> have shown skilled techniques in designing self-supporting filters, and the self-supporting sensitivity of each element was expressed in terms of the densities information of the layers below it. The approach has its great freedom in computing self-supporting structure, but also requires a further improvement to overcome it inhabitation of parallel processing, as also explained by the authors <cit.>.
    
    Derivations of the sensitivities of the objective function(̧)or of the volume constraintV(), involved in (<ref>), are totally the same as those done in previous studies <cit.>. Integrating these sensitivities with thickness control can be achieved using the Heaviside filter <cit.>, as will be demonstrated in Section <ref>. Details are not further explained here.
    
    
    
     §.§ Discrete convolution for efficient unsupported element detection
    
    Computing the sensitivity (<ref>) of the self-supporting constraints requires detecting the set_Uof all the unsupported elements. They can be easily detected via enumerating all the discrete elementseofnot satisfying the property in (<ref>) or (<ref>). Such process is however very time-consuming. In order to further accelerate this process, a novel convolution operator is designed for such detections, which acceleration the process of detecting the unsupported element with a speedup of 100 times as compared with the approach of directly enumerating them element-wise. 
    
    pGiven a discrete structureof sizeN× Min 2D, we can see from (<ref>) that an elemente(n,m)∈is supported if the summation of the densities of its supporting elements is larger than zero, or specifically,
    
        ∑_n-1≤ r≤ n+1ρ(r,m-1)> 0.
