Topology optimization incorporating a passageway for powder removal in designs for additive manufacturing

In powder-based additive manufacturing, the unused powder must be removed after printing. Topology optimization has been applied to designs for additive manufacturing, which may lead to designs with enclosed voids, where the powder will be trapped inside during printing. A topology optimization method incorporating a powder removal passageway is developed to avoid the powder being trapped inside the structure. The passageway is generated by connecting the entrance, all voids, and the exit sequentially. Each void is limited to have only one pair of inlet and outlet to guarantee a single-path flow to facilitate powder removal after the additive manufacturing. The path of the passageway is optimized to minimize its influence on structural stiffness. The proposed optimization method was applied to two practical case studies where the powder removal passageways were generated successfully.


Introduction
Additive manufacturing has experienced rapid development since the creation of the first 3D printer (Hull 1986). There have been various additive manufacturing techniques, such as fused deposition modeling (FDM), stereolithography (SLA), selective laser melting (SLM), and selective laser sintering (SLS). The principle of most additive manufacturing methods is building a component in successive layers, which effectively differentiates additive manufacturing from traditional subtractive manufacturing methods. The unique working principle of additive manufacturing has made it suitable for building components with complex structures (Gao et al. 2015), such as truss lattice structure (Rosen 2007;Wang et al. 2018;Yu et al. 2020), irregular shapes (Thompson et al. 2016), and internal features (Thompson et al. 2015;Liu et al. 2020). For powder-based printing methods such as SLM and SLS, the printing is conducted in a chamber filled with powder, and the unused powder needs to be removed after printing (Kruth et al. 2004(Kruth et al. , 2005Yap et al. 2015).
There are a few cleaning methods, including manual brushing, vibrations, plasma shockwaves cleaning (Dunbar et al. 2007), ultrasound cleaning (Mittal and Jaiswal 2015), and compressed air flushing (Vayre et al. 2013). For components with internal features, an appropriate design should be used so that the powder can be effectively removed from internal voids.
Design for additive manufacturing is becoming an important step in additive manufacturing (Gao et al. 2015). Topology optimization has played an important role during the design process because of its capability of reducing component weight while maintaining the structural performance as much as possible . The designs from topology optimization usually have irregular and complex shapes, which makes additive manufacturing a suitable choice for building the component (Meng et al. 2020). There have been various topology optimization methods, such as homogenization method (Bendsøe and Kikuchi 1988), evolutionary structural optimization (ESO) method (Xie and Steven 1993, 1 3 232 Page 2 of 17 1997; Huang and Xie 2010a), bi-directional evolutionary structural optimization (BESO) (Querin et al. 1998;Yang et al. 1999;Huang et al. , 2014Huang and Xie 2010b;Zuo and Xie 2015), solid isotropic material with penalization (SIMP) (Bendsøe 1989;Zhou and Rozvany 1991;Rozvany et al. 1992), level-set method (Osher and Sethian 1988;Wang et al. 2003;Allaire et al. 2004), and element volume fraction-based methods (Guo 2014;Da et al. 2018;Fu et al. 2020b;Huang 2020).
Additive manufacturing has removed many limitations of traditional manufacturing, but some remain (Liu and Ma 2016). One of the constraints is the limit on the overhang length (Calignano 2014). Additive manufacturing can build complex overhanging structures, while the overhang should be less than a certain threshold to avoid the collapse of the overhanging section. There is substantial research attempting to resolve this issue. For FDM, the sloping wall structure ) and automatically generated stable scaffolding (Dumas et al. 2014) can be used to support the overhangs. For SLM, an "IY"-shaped support was developed to reduce the contact area (Gan and Wong 2016), and optimization was applied to the support structures to reduce the residual stress. There are also some efforts made to reduce the support structures during topology optimization. This includes the modification of large voids into multiple support-free voids (Leary et al. 2014), application of undercut and overhang angle control in the topology optimization to reduce support structures (Qian 2017), application of void size control in the topology optimization to obtain supportfree voids (Stojanov et al. 2017), and application of additive manufacturing filter in the topology optimization method to generate self-supporting designs (Fu et al. 2019(Fu et al. , 2020a. The support structures are removed by these methods, while more voids are created with the reduction of void size. Another constraint specifically for powder-based additive manufacturing is the requirement for powder removal. A feasible design is needed to ensure the powder can be removed after printing. There have been some studies focusing on the complexity control in topology optimization, which is usually characterized by the distribution and geometries of voids. Zhang et al. (2017) controlled the structural complexity through the constraints on multiple level-set functions in a level-set-based topology optimization. Zhao et al. (2020) applied void number and size constraint in evolutionary topology optimization. They used algorithms from graph theory and set theory to describe and adjust the status of voids formed by elements. There were also some studies directly working on the no-enclosed-voids constraint in topology optimization for additive manufacturing. Liu et al. (2015) developed a virtual temperature method, which served as a simply connected constraint in topology optimization. The voids were assigned with a virtual heating material, and the solid areas were assigned with a virtual insulating material so that the structural connectivity was related to the system temperature. By applying a maximum temperature constraint, the enclosed voids were pushed towards the boundaries and finally eliminated. This method was further developed to a virtual scalar field method for the connectivity constraint in topology optimization (Li et al. 2016), and applied to the optimization of cast parts . Hu et al. (2017) applied the virtual temperature method in the optimization of a space telescope to eliminate enclosed voids in the final design. Luo et al. (2020) further proposed a non-linear virtual temperature method which limited the temperature to a specified range in the structure. The virtual temperature methods eliminated the enclosed voids by pushing them towards the boundary, while the geometries of the optimized designs could be substantially changed, and it might be hard to handle the models where a long pathway is needed. The application of the virtual temperature method also requires further calculations of the virtual scalar field. Gaynor and Johnson (2020) combined the overhang projection method and void projection method to eliminate the occluded voids in the designs. Their approach was claimed to be able to generate long exit pathways while not able to allow for certain serpentine void pathways. Zhou and Zhang (2019) used voids as design features and limited all void centers outside the design domain using side constraint method in a level-set-based topology optimization, while the optimized designs could largely differ from the designs without constraints. Xiong et al. (2020) applied connectivity control in the evolutionary structural optimization. The enclosed voids are eliminated by generating tunnels between voids and boundary, while the flow path for the powder removal is not considered.
In summary, there have been some work trying to eliminate the enclosed voids in the topology optimization process. They are mostly applied to geometries with the freedom to alter the outside shape of the original geometry. Little work has been conducted on re-designing existing solid components under the constraint of retaining the outside surface of the original design. Furthermore, little work has been done on generating designs with the consideration of extracting remaining powder during the powder-based additive manufacturing. In order to maintain the external surface geometries while generating internal voids and, therefore, to provide a path for powder removal after the printing, a method incorporating a powder removal passageway to link enclosed voids is proposed. During the topology optimization, the structural connectivity is evaluated by identifying the enclosed voids in the structure, and then a powder removal passageway is generated to connect the enclosed voids to the external boundary. The passageway will cross through the entrance, all voids, and then the exit one by one to guarantee a continuous path for powder removal. The potential regions for entrance and exit can be assigned to avoid possible intrusion into functional features. The path is also optimized to reduce the increase in the objective compliance of the structure.

Methodology
In order to eliminate enclosed voids in the designs for additive manufacturing to facilitate powder removal, a powder removal passageway is added to the topology optimization method. The workflow of the whole method is shown in Fig. 1. The optimization starts with the initialization of the model and parameters, which are read and stored into a database to facilitate further usage. The following step is to carry out finite element analysis (FEA) and calculate sensitivities. The elements are ranked based on their sensitivities. The ones with sensitivities higher than the threshold are treated as solid elements, and the ones with sensitivities lower than the threshold are deactivated and treated as void elements.
This applies to the whole model including the passageway generated in the last iteration. Internal voids, formed by clusters of void elements, are then identified based on the connectivity status between void elements through depthfirst search (DFS). A connection graph is built based on the connectivity among elements and sensitivities of elements. The minimum cost path, which refers to the elements path with minimum overall sensitivities, between each pair of voids is obtained using Dijkstra's algorithm. The minimum cost path through all the voids is then obtained by a genetic algorithm. Then the passageway is generated and the whole design is updated. The optimization stops when the halting condition is met.

Initializing the optimization
To prepare for the topology optimization method, the model is discretized into finite elements. A database is built by scanning the model and extracting the data. A model can Fig. 1 Flowchart of the topology optimization process incorporating a passageway for powder removal contain several parts including the design regions and nondesign regions. The setup of multiple parts in a model facilitates the optimization of practical models where there can be several components inside. The part can contain many elements, the neighboring status of each element is checked and recorded in the database to facilitate the passageway generation.

Calculating sensitivity values
The BESO method (Huang and Xie 2010b) is the topology optimization method used for the present optimization task. For a given model with a specified design space, there will be many finite elements. The existence status of element i is represented by x i , where x i = 0 stands for void and x i = 1 stands for solid. The method uses a hard-kill version of BESO to facilitate the void identification. It can also be switched to a soft-kill version of BESO. For a static, linear elastic structure, the optimization task to minimize compliance is written as follows: where C is the objective function, f is the force vector, u is the displacement vector, V* is the prescribed volume constraint, V i is the volume of element i, and N is the total number of elements. The sensitivity of element i is calculated as (Huang and Xie 2010a) follows: where 0 i stands for the stiffness matrix of element i, i is the displacement vector of element i. For different x i value, the sensitivity number becomes It should be noted that Eqs. (2) and (3) hold the same for models with body force (Huang and Xie 2010b). The (1) sensitivity values are then smoothed to obtain a meshindependent solution (Díaz and Sigmund 1995;Sigmund and Petersson 1998;Sigmund 2007). The nodal sensitivity numbers are firstly calculated by averaging the elemental sensitivity numbers: where M is the number of elements attached to the jth node, w i is the weight factor of the ith element and where r ij is the distance between the center of the ith element and jth node. It shows that closer element has larger effect on the nodal sensitivity number. The nodal sensitivity numbers will then be converted into smoothed elemental sensitivity numbers by a filter with radius r min . If a specified passageway radius is used, the filter radius should be larger than the passageway radius to guarantee effective update of the passageway. Nodes within the distance of r min to the center of one element will be used to calculate the smoothed sensitivity of this element: where K is the total number of nodes contributing to the smoothing, w r ij is a weight factor defined by The sensitivity value is further averaged with its historical value through iterations to stabilize the convergence of the objective function .
where k is the current iteration, and k − 1 is the previous iteration. The processed sensitivities are used for ranking the elements. Each element's status is modified based on the rank of the element sensitivity and a threshold number decided by the target volume for the next iteration V k+1 . In each iteration, the target volume is where ER is the evolution ratio that represents the percentage of volume variation in each iteration. When the volume constraint is satisfied, the volume will be kept constant for the remaining iterations as V k+1 = V * . The elements with sensitivities higher than the threshold are treated as solid and the elements ranked below the threshold are treated as void.
The convergence condition from Huang and Xie (2009) is used: where k stands for the current iteration, N is an integer representing 2 N iterations used for the convergence check, is allowable error for convergence. A general recommendation is N equals to 5 and equals to 0.01% (Huang and Xie 2010b), while in the case studies, a few more iterations after the convergence are included to facilitate the comparison of results. For models with uniform cross-section requirements, a grouping operation method developed by Stojanov et al. (2016Stojanov et al. ( , 2017 is applied to make sure the BESO design has a uniform cross-section. The grouping operation method can be applied to a quasi-2D model with a constant cross-section or a revolved model with axisymmetric constraint. The elements along the thickness or circumferential direction are treated as one group, as shown in Fig. 2. The elements perpendicular to the grouping direction are treated as the cross-section, which should be specified manually before the optimization. The sensitivity of the group is calculated by averaging through the elements inside. The ranking, deactivating, and activating operations are applied to the whole group rather than a single element. This can guarantee a design with a uniform cross-section or axisymmetric geometry regardless of the loading condition. Illustrations of a discretized model a with solid elements (orange) and void elements (white), b with enclosed voids (green), potential region for entrance (red) and potential region for exit (blue), c the internal void elements and boundary elements (gray) assigned with high sensitivities 232 Page 6 of 17

Creating a passageway through all voids
The finite element analysis, sensitivities calculation, sensitivities filtering, elements grouping, and elements ranking processes detailed in the previous sections are applied to the model. The existence status of the elements is updated based on the ranking results, and the elements with lower sensitivities are deactivated. These deactivated elements form internal voids, as shown in a sample model in Fig. 3. The passageway building module assures the connectivity between the enclosed voids and external space. This includes identifying voids, searching for optimal paths between void pairs, and designing an optimal path crossing through all voids, as presented in the following sections.

Identifying voids
The initial design from previous step contains both solid elements and void elements. The enclosed voids formed by the void elements are to be connected to the external space.
To create a passageway through all voids, the enclosed voids are firstly identified by dividing void elements into groups.
The whole model consisting of finite elements is then treated as a graph G = (V, E), where the elements are represented by vertices V and the connection status between elements are represented by edges E. A void can be viewed as a subgraph S = V S , E S , where for each x, y ∈ V S , there is a path from x to y, and vice-versa. The task of identifying voids in the model is equivalent to finding a strongly connected component in the graph. An algorithm from Pearce (Pearce 2005) is applied to find the strongly connected void elements. The basic idea of the algorithm is using a depth-first search (DFS) algorithm to transverse through the graph and dynamically record the status of the visited vertices. The first step of void identification is picking out all void elements by checking the existence status of elements. These void elements are passed to the following searching steps.
After that, the void is identified by searching through the connected elements. In a 3D discretized model, two elements are viewed as connected if they share a common face. The connectivity status between elements is updated and stored iteratively in the database, as mentioned in Sect. 2.1. When the searching algorithm is applied, the connectivity status can be directly read from the database. A sample void with five void elements is shown in Fig. 4. The DFS algorithm will firstly search through the depth, visiting elements 1, 2, 3, and 4, then go back to element 5 in the side branch. All the visited elements are labeled as visited and recorded in a cluster, which represents a void. After that, the DFS algorithm is switched to the unvisited elements to find another void. The searching stops when all void elements are labeled as visited.
A boundary element set is defined before the optimization to indicate the boundary of the model. It is grouped as an element item set and stored in the database. The void elements connected to the boundary set is treated as belonging to external space, which should not be viewed as internal void elements. All void elements clusters are checked against the boundary. If any of the elements in a cluster is connected to or belong to the boundary, the whole cluster will be removed and excluded from the following operations. The remaining clusters represent separated enclosed voids consisting of void elements. The enclosed voids of the sample model are highlighted in green, as shown in Fig. 3b. Potential regions for entrance and exit can be artificially set to the chosen locations, as highlighted in red and blue in Fig. 3b. This would provide a choice to the users to avoid possible influence on certain functionality features in the designs. For simple models without such requirements, the entrance region and exit region can be set to the whole external boundary.

Optimizing paths between void pairs
The whole model consisting of finite elements is viewed as a graph G = (V, E). The voids are treated as the subgraphs that have been identified in the previous section. The passageway can be viewed as a path P = [v 0, v 1, …, v x ] connecting and crossing through the subgraphs one by one. The search for the optimal path is divided into two steps. The first step is finding out the optimal path between each pair of two voids, which is detailed in the present section. The second step is finding out the optimal sequence for the path crossing through the voids. There can be many possible paths P between any pair of two voids. A path with a minimum cost of structural stiffness is expected from these possibilities. The optimal paths between void pairs are prepared as the potential options for the following step. Dijkstra's algorithm is then applied to find the minimum cost paths.
Firstly, an adjacency matrix is built across all elements. The connectivity status of elements and the costs of these Fig. 4 The DFS algorithm searches through the void elements in a sample void connections are stored in the matrix. The connection cost between two neighboring elements is calculated as the average of their sensitivities: where e v i v j is viewed as the length of the edge connecting vertex i and j in the graph, i and j are the sensitivities of corresponding elements i and j, respectively. This converts the minimum cost path problem to the shortest path problem. The objective function is optimized by searching for the shortest path.
The passageway size can influence the performance of powder removal. An experimental study by Hunter et al. (2020) shows that in average around 63% powder is removed from their designed geometry with 1-mm-diameter passageway, 78% powder is removed from the same geometry with 2-mm-diameter passageway, and 75% powder is removed from the same geometry with 3-mm-diameter passageway. Both ultrasonic polish + ultrasonic bath and ultrasonic polish + vacuum boil were used to obtain the data. It shows the importance of providing a choice for the passageway size. For the situation where a specified passageway radius is needed, the element sensitivities are firstly averaged through all elements within the distance of a given radius using Eq. (10) when building the adjacency graph.
The internal void elements and all boundary elements are then assigned with high sensitivities, as shown in Fig. 3c. In topology optimization, the internal void elements were deactivated because they had lower sensitivities. During the search of the low-cost path, these elements have higher possibilities of being selected. It means a path connecting two voids could possibly cross through a third void, which will lead to multiple inlets and outlets in a single void. This could influence the flowability of the passageway and the effectiveness of the powder removal. To avoid such a scenario, the sensitivities of void elements are replaced by a large value α up , which is ten times the maximum sensitivity in the model. The path between two voids will tend to bypass the other voids since the void elements are assigned with these higher sensitivities. Similarly, the boundary elements are also assigned with a large sensitivity number to avoid the path connecting to the boundary unexpectedly. The path between enclosed voids and external space should only cross through the specified entrance and exit. It should be noted that the adjustment of sensitivities only occurs within the adjacency matrix for path searching, and the elements ranking for topology optimization is not influenced. Dijkstra's algorithm (1959) is applied to find the minimum cost path between two voids. The pseudocode of the algorithm is shown in the Appendix A. The algorithm starts from the source vertex and checks the distances of all its neighbors to the source. After that, it proceeds to the closest vertex and again checks all neighbors of this vertex. The distance of each vertex to the source and the path between them is dynamically updated. When the graph is traversed, the shortest path from the source to each vertex is obtained. The algorithm is applied to all void elements inside the enclosed voids. The connections between each pair of two voids are gathered, and the one with the lowest cost is the target path. It should be noted that the Dijkstra's algorithm only works with positive weights. For optimization with negative sensitivities, the values can be adjusted by data processing operations, such as adding a constant to all sensitivities when preparing the graph to avoid negative weights. This operation is applied to the graph only and has no effect on other parts of the optimization.
The search is also applied to the entrance and exit set of the model, which are specified before the optimization. The minimum cost path between the entrance/exit and all voids are also obtained. The entrance and exit can be constrained to appropriate positions. This effectively assures that the passageway openings will not influence the functionality of critical surfaces.

Optimizing sequence for the passageway crossing through all voids
The previous section delivers the optimal paths between any two voids. Then an optimal sequence for the passageway Fig. 5 The map of costs between the entrance, exit, and all voids in the sample model crossing through the voids is required. The passageway should start from the entrance, cross through the voids sequentially, then end up at the exit. A sample model is illustrated in Fig. 3b, with the entrance highlighted in red, voids highlighted in green, and exit highlighted in blue. The costs between every two items are obtained from the previous modules, and the cost map of a sample model is shown in Fig. 5. An optimal passageway should have a minimum cost and visit each void only once. This is equivalent to finding the shortest Hamiltonian path in the costs map. It is well known that this problem is NP-complete. Therefore, a genetic algorithm (Taylor 1994) is applied to solve this problem, which is briefly presented in Appendix B.
The processes in previous sections are repeated until the stopping condition of the optimization is met. The final design from the optimization will contain voids created by topology optimization, and one passageway that crosses through all the voids. The passageway is updated when the voids are updated. It is evolved along with the topology optimization progress to achieve an optimized path. Post-processing methods are applied to the voids and the passageway to obtain a final design with smooth boundaries. It is worth mentioning that although this proposed method can ensure a passageway crossing through all the voids, it cannot completely avoid dead-corners in the voids, where the powder may reside and therefore may be hard to be removed from. The choice of entrance and exit locations may mitigate such situations and removing powder from dead-corners may also be resolved by using novel powder removing techniques. To completely avoid dead-corners including how to unambiguously define a dead-corner requires further study.

A jet engine impeller
A jet engine impeller is used as a case to demonstrate the topology optimization with a powder removal passageway. The impeller, which is a critical component in the engine compressor, is used to compress the air through fast rotation.
In the aerospace industry, it is well known that the component's weight can substantially influence the fuel cost in the operation of aircraft, which is especially true for a rotating component. In this case study, topology optimization is used to reduce the impeller's weight, and a powder removal passageway is included to make sure the powder in the internal voids can be removed after printing.
The original model of the impeller is shown in Fig. 6a, where centrifugal loads are applied to the whole model, pointing outside along the radial direction. 1/10 of the hub, as shown in Fig. 6b, is used for optimization since an axisymmetric constraint is applied to the model. The inner hub of the impeller is set as the design region, and the fins are not touched. The potential regions for entrance and exit of the passageway are both set at the inner side of the impeller, as highlighted in Fig. 6b. This is because the outer surface is designed for compressing air, and the passageway outlet over such surface will influence the performance of the compressor.
Topology optimization is performed on the section of the impeller. The section model has a volume of around 2.9 × 10 4 mm 3 , and 0.5 mm cubic elements are used in the model. A 95% volume ratio of the design region is used as the optimization constraint to obtain a design with many small voids so as to showcase the capability of generating passageway through multiple voids of the method. The evolution ratio is set to 0.3% and the filter radius is set to 1 mm. The optimization design without connectivity constraint is shown in Fig. 7a. The voids are enclosed inside the structure, where the powder cannot be removed. The optimization design with the powder removal passageway is shown in Fig. 7b. It is noted that the voids in Fig. 7b are slightly different from those in Fig. 7a due to the interaction between the passageway and voids in the proposed optimization. The passageway starts from the entrance, crosses through the voids, and then ends at the exit. It should be noted that the passageway is not revolved around the circumference, but only patterned along with the 1/10 slice in this case study. The optimization histories of both models are shown in Fig. 8, where only 0.5% of the overall stiffness is lost from the optimization. The compliance of the optimized impeller with passageway is 100.004% of the optimized impeller without passageway. The section design is smoothed and converted back to the full model by circular patterning of the single section. The smoothed impeller with the powder removal passageway is shown in Fig. 9. The passageway is patterned ten times along with the section to obtain the full model.

A rocket fin
A rocket fin is the second case study for demonstrating the optimization with a powder removal passageway. A fin is used at the tail of a rocket to improve the stability of the flight attitude of the rocket. Weight reduction is crucial to the performance of a rocket. Topology optimization with a powder removal passageway is applied to provide a lightweight design with no enclosed voids. The original fin is illustrated in Fig. 10. A traction load is applied on both side surfaces of the fin, pointing towards the tail horizontally to simulate the traction from airflow.
The side surfaces of the rocket fin are used to guide the airflow, where no modification should be made. Therefore, these surfaces are set as the non-design region in the topology optimization, and only the internal materials are used as the design region. The entrance and exit are set to the top and bottom surface of the fin, as shown in Fig. 10a, where the passageway openings are less likely to influence the functionality of the fin. The volume of the whole model is around 1.4 × 10 5 mm 3 , and 2 mm elements are applied to the model. A 70% volume ratio is used for optimizations without and with powder removal passageway. The evolution ratio is set to 0.5% and the filter radius is set to 3 mm. The topology optimization design without connectivity constraint is shown in Fig. 11b. The fin model is cut through the middle plane to show the internal structures. The topology optimization design with powder removal passageway is shown in Fig. 11c. The passageway starts from the entrance at the top surface, crosses through the two voids, then ends at the exit at the bottom surface. Similarly, the voids in Fig. 11c are different from those in Fig. 11b due to the interaction between the passageway and voids during the optimization. It is worth noting that an alternative solution to the case is having two separate passageways connecting each of the two voids, although this is out of the scope of the current study. It would be an interesting future development of the method. The optimization histories are shown in Fig. 12, where approx. 3.5% stiffness is lost from the deletion of elements. It is worth noting that there is a sudden jump in the compliance at the beginning of the optimization with a passageway, which is caused by the appearance of a long passageway connecting the entrance and the exit. The compliance of the optimized fin with passageway is 100.095% of the fin without passageway. The final design, where the voids and passageways are post-processed by smoothing filter and spline fitting, is shown in Fig. 12d.

A rocket fin with a specified passageway size
The proposed method provides the choice of passageway size. In previous case studies, the passageway is generated with only one element width. In this case study, the same rocket fin used in Sect. 3.2 is optimized with a specified passageway radius. The average size of elements used in the design region of the rocket fin is around 2 mm. The evolution ratio is set to 0.5% and the filter radius is set to 3 mm. The passageway radius is set as 2.25 mm to cover the larger elements. For each element, the sensitivity is averaged through the elements within the region of 2.25 mm radius. All the other parameters including loadings, boundary conditions, target volume ratio, and are set to be same as the fin model in Sect. 3.2. The number of iterations is artificially set to 80 to facilitate the visualization of optimization history and comparison of computational time.  The topology optimized design is shown in Fig. 13a. The elements within the distance of 2.25 mm to the elements along the selected path are removed to generate the final passageway. The elements in the model have various size, therefore irregular shapes are generated in the passageway. Some staircases-like features are algo generated from the passageway radius filter. The optimized model is then smoothed, as shown in Fig. 13b. The passageway has a much larger size compared to the model in Fig. 11d. The optimization history of the fin model is shown in Fig. 14. There is a big jump in the normalized compliance at iteration 1. This is because the passageway is connecting the entrance and exit through only a small void. The long and wide passageway leads to a large loss of compliance. The proposed topology optimization method with passageway generation function spends around 13% more time compared to a topology optimization without a passageway in this case.
To demonstrate the printing feasibility of the optimized designs based on the proposed method, an impeller and a fin were printed in steel by using selective laser melting based on preliminary designs, which are presented in Appendix C.

Conclusions
A topology optimization method incorporating an internal passageway for powder removal is presented in this manuscript to deliver feasible designs for powder-based additive manufacturing. The passageway starts from the entrance, crosses through all enclosed voids, and ends at the exit to provide a single-pass path for compressed air to flush out the unused powder. The shortest path algorithm and genetic algorithm are applied to find the path with minimal influence on the structural stiffness. The potential regions for entrance and exit can be artificially assigned on the model before the optimization, which is different from the existing research where the locations of outlets are not controlled. This can be useful, especially for practical components where the design space should be carefully selected to avoid compromising the functionality of design features, as presented in the two cases. Each void is limited to only one pair of inlet and outlet to guarantee the air flowability through the voids. Two practical designs, an impeller for a jet engine and a fin for a rocket tail, were used to demonstrate the method. Designs with internal voids and a powder removal passageway were successfully generated using the presented method in the two case studies.

Appendix B
The workflow of a typical genetic algorithm is shown in Fig. 15. The algorithm starts by creating an initial population of chromosomes, which represents the sequence of the passageway crossing through the voids. Assuming there are six voids, the initial population of the chromosomes could be the ones shown in the figure. These chromosomes are generated randomly. During the random selection, the entrance is always at the beginning, and the exit is always at the end. The fitness of the chromosome is defined by the overall cost of the sequence. This is obtained by summing up the sensitivities of the elements along the whole passageway. After the evaluation and ranking of the fitness for all chromosomes, only the ones within the lowest x% are kept in the population. Then the crossover and mutation are applied to these chromosomes to create new chromosomes to form a new population. The crossover is performed by combining the sequence of two chromosomes. For the two chromosomes shown in Fig. B, the initial sequences are Ena-b-c-d-e-f-Ex and En-d-f-c-e-a-b-Ex. We generate a random slice of index 2-3-4 from the whole sequence of index 1-2-3-4-5-6. The corresponding slices b-c-d from the first chromosome and f-c-e from the second chromosome are firstly switched. The remaining genes are assigned to the remaining positions following the original sequence. For example, after excluding the genes f, c, and e from the first chromosome, the remaining genes are a, b, and d. These genes are assigned to position 1, 5, and 6, respectively. Therefore, the first chromosome becomes En-a-f-c-e-b-d-Ex, and the second chromosome becomes En-f-b-c-d-e-a-Ex after the crossover. The mutation is performed by switching the position of two voids randomly. For example, chromosome En-a-b-c-d-e-f-Ex could be changed to En-a-e-c-d-bf-Ex. This process is repeated until the fitness converges or the maximum number of iterations is met.