Fluid topology optimization and additive manufacturing of a liquid atomizer using an extensive number of grid points

Additive manufacturing (AM) can fabricate complicated shapes and is useful when manufacturing topology optimized shapes. Fluid parts often consists of three dimensional curves that are suitable for AM fabrication. However, the application of fluid topology optimization and AM has not been investigated yet. However, modeling and solving an optimization problem have not been investigated for a real industrial fluid topology optimization problem of AM parts with tiny channels, i.e., a liquid atomizer which is equipped with an aero-engine fuel injector. In order to reduce computation time, which is an important issue in real industrial problem, the instantaneous sensitivity approximation method is used as a topology optimization method. The optimized part exhibited a reduction in pressure loss compared to that of a conventional part.


Introduction
Additive manufacturing (AM) techniques can be used to overcome the limitations of conventional manufacturing approaches and fabricate complex shapes. In general, conventional machining tools enforce several restrictions on the manufacturing of complex parts. Specifically, internal channels are difficult to be manufactured because machining tools can drill through only one direction. In the context of hydrodynamics, it is desirable to reduce the pressure loss in fluid applications by using smoothly curved channels instead of right-angled channels. However, threedimensional curved channels cannot be easily manufactured using conventional tools.
AM is also useful when fabricating shapes that are obtained by topology optimization methods which usually consist of free smooth curves [12,22,25,44]. Many research that couples AM and fluid topology optimization has Kazuo Yonekura yonekura@struct.t.u-tokyo.ac.jp 1 The University of Tokyo, 7-3-1, Hongo, Bunkyo-ku, Tokyo, 113-8656 Japan 2 IHI Corporation, 1, Shin-nakahara-cho, Isogo-ku, Yokohama, Kanagawa 235-8501 Japan been conducted [21,42]. Jahan et al. [15] used topology optimization for conformal cooling channels that is fabricated by a plastic 3D printer. Ghosh et al. [14] optimized and manufactured a U-bend channel, and experimentally validated the pressure loss. However, real industrial applications have not been well reported. Especially, a topology optimization of tiny channels that needs large computation time has not been investigated. One of difficulties of modeling and solving such real application is that the fluid computation needs extensive number of grid points which leads to a long computation time. In this study, we focus on a real industrial application, i.e., a fuel atomizer. The fluid loss is minimized using the fluid topology optimization method called instantaneous sensitivity approximation (iSA) method, whose computation cost is small, and fabricate it by metal AM machine to discuss its usefulness. AM techniques can be used to fabricate threedimensional curved channels, for instance, in the case of fuel nozzles. General Electric Company re-designed and fabricated a fuel nozzle using AM, and this nozzle has been employed in an operational aero-engine [10]. A fuel nozzle is a fluid device used to feed fuel into a combustor. As mentioned previously, if the internal channels in a nozzle are fabricated using machining tools, the flow channels consist of narrow straight holes and involve several rightangled corners. Consequently, a large fluid loss occurs, and a large actuator is necessary, which increase the weight of aircraft engines. It is studied in literature if AM is applicable for liquid atomizers. Manufacture [1] studied that the characterization of an atomizer manufactured by metal AM, and reported that the atomizer successfully made a uniform spray pattern. The computational fluid dynamics (CFD) method in applications to aero engines are studied in literature [9,11,32], and is used in a design process.
To maximize the advantages of AM, it is desirable to optimize the topology of the parts [6]. Several researchers have reported on the manufacturing of topology-optimized parts considering structural aspects such as the compliance [21,24,30,41]. Nevertheless, AM involves certain manufacturing limitations. For instance, manufacturing cannot be performed in a closed room as the metal powder in the room cannot be removed. Moreover, the angle of the overhang area must be less than a certain threshold; in other words, if the overhang angle is too steep, the part cannot be manufactured using AM, and a support structure is necessary to support the overhang area. Several researchers have attempted to overcome these limitations pertaining to the powder removal [30] and the overhang angle and support structure [13,18,23,33]. The necessity of support materials is especially important in fluid channels because internal support structures cannot be removed after builds. Behrou et al. [5] optimized fluid parts without internal supports, but the optimized part was not actually fabricated. Furthermore, the heat exchange problem has been considered to optimize lattice cooling devices [31], and some researchers reported on the topology optimization for nano-photonics as well [16]. Behrou et al. [5] proposed a topology optimization method for a fluidic component that is free of internal support. However, in those literature, textbook examples are studied, and industrial application of such fluid topology optimization has not been studied. Especially, manufacturing the optimized shape by metal AM has not been reported. One of the difficulties in applying AM and topology optimization to an industrial problem is that the industrial parts are complicated and hence need a massive number of grids for computation. In the present study, the fluid topology optimization method is applied to the industrial application, which is the liquid atomizer, where a massive number of grids are necessary to compute fluid flow. Moreover, the optimal shape is fabricated by metal AM.
To optimize shapes while maintaining the fluid performance, several topology optimization methods for fluids have been reported [3,7,20,28,35,36]. However, a limitation of fluid topology optimization is the large computation time required. In general, industrial applications involve several degrees of freedom, and thus, the computation cost must be reduced. Structural topology optimization has been applied to large-scale industrial problems. The state-of-the-art is 1.1 billion voxel in [2]. In fluid topology optimization, the lattice Boltzmann method (LBM) based optimization method is expected to be comparatively time efficient. Using this approach, [34] solved a problem involving 50 million grids. A topology optimization method using temporal sensitivity, known as the instantaneous sensitivity approximation (iSA) method, was developed to enable fast computation [37][38][39][40]. It is reported in [37] that the computational times of iSA to solve around 3000-grid-pointproblems were less than 11 s. The conventional optimization approach would take more than few minutes to solve such problems. Hence, the iSA method is drastically faster than conventional methods.
Considering these aspects, in this study, we optimized the internal channel path of a liquid atomizer to reduce the pressure loss compared to that in a conventional part manufactured via machining. In order to simulate the flow in the channel, about half a billion grids were necessary for this problem, which is extensively large. To conduct optimization in a reasonable time, the iSA approach was employed. Furthermore, the optimal shape was fabricated using AM approaches to demonstrate the manufacturing feasibility. The fabricated parts were polished to ensure surface roughness. A schematic diagram of design and manufacturing the part is described in Fig. 1.
The contributions of the present paper are reporting fluid topology optimization method to an real industrial problem, overcome the issue of computation time by iSA method, and reporting polishment after fabrication by AM machine.
The remaining paper is organized as follows: Section 2 describes the topology optimization method and the liquid atomizer optimization. Section 3 describes the numerical fluid analysis performed to verify the optimization result and the fabrication of the liquid atomizer. The concluding remarks are presented in Section 4.

Liquid atomizer
Gus turbine engines compress air, burn fuel, and extract energy by using turbines. Fuel nozzles are used in gas turbine engines to atomize fuel and inject it into a combustor to realize fuel combustion under favorable conditions. The mechanism of liquid atomization has been extensively investigated [4,8,19]. In this study, the objective was to optimize the swirl atomizer fuel channel [4], wherein liquids flows through the wall inside a cylinder. In general, such atomizer are small, with the diameter of the cylinder being typically less than 1 in. The liquid channel before the swirl chamber has a complex shape and is usually manufactured via machining. Specifically, the atomizer is separated into several pieces and assembled into one part. The shape of the internal channel is thus limited owing to the machining process, leading to a large fluid loss. However, using AM, more complex shapes can be manufactured. Consequently, it is desirable to optimize the flow channel.
We set the design domain as shown in Fig. 2, and the design problem was formulated as follows: The fuel is provided from the bottom side through the duct and input into the pool. Subsequently, the fuel flows through four channels into the swirl chamber. The swirl chamber is appropriately designed considering the atomization characteristics, and its dimensions cannot be modified. The number of channels is fixed as well. Atomization process needs a large energy dissipation. If the whole domain including channels and swirl pool is set as the design domain, the minimization of energy dissipation in the channels and the best atomization characteristics are incompatible. Therefore, only the shape of the four channels in the design domain was optimized while minimizing the energy loss.

Optimization method
In the case of such an extremely large number of grids, the computation time is unreasonably large. Therefore, we used the iSA method [37][38][39] that can enable a reasonably fast computation. Specifically, we used the level-set function coupled iSA method [39] to obtain a clear representation of the fluid and solid regions. The program was implemented in C++, and openMPI was used for parallel computation. In the iSA method, the flow field is generally computed using LBM [29], which is a numerical method for computational fluid dynamics and is widely applied in various fields including AM metal melting analysis [17,26].
The flow chart of the iSA method is described in Fig. 3. The method updates flow field and design variables simultaneously. The flow field solved by LBM is unsteady, and the sensitivity to update the design variables is computed using the unsteady flow information. Hence, the sensitivity is called an instantaneous sensitivity.
The optimization objective was to reduce the fluid energy loss, and the optimization problem was formulated as follows: where u and φ represent the flow velocity and a design variable, respectively. 0 ≤ φ ≤ 1 is defined on each computation grid. The function g is defined as g(y) = 1 (y ≥ 0.5), 0 (y < 0.5). If φ < 0.5, the grid pertains to a solid region; otherwise, the grid pertains to a fluid region. Hence, a part is represented by a voxel-like shape. The objective function pertains to the energy dissipation of the flow Eq. 1b inside the fluid domain Ω f , where ν and ε represent the viscosity and strain velocity tensor, respectively. As indicated in Eq. eq.opt0.e1, the velocity is computed based on the LBM. The volume of the fluid region is restricted, as given in Eq. 1d. Although this volume constraint is not necessary from a practical point of view, it is vital for the optimization formulation. We setV = 0.3Ω, where Ω is the volume of the overall design domain. The optimization process in most fluid topology optimization methods involves the following steps: First, φ is fixed, and the flow field u is computed by solving the fluid dynamics equations. Next, the sensitivity for φ is computed, and φ is updated. The flow field is determined again considering the new φ, and these steps are repeated until the optimization is complete. In contrast, in the iSA method, the flow field variables and design variables are updated simultaneously. The algorithm solves the unsteady flow equations, i.e., the lattice Boltzmann equations, for only one time step. Subsequently, the sensitivity is determined using the flow information pertaining to the yet non-steady state. This sensitivity is termed as the instantaneous sensitivity. φ is updated using the instantaneous sensitivity. Subsequently, the flow equations are evolved for another time step. The algorithm terminates when the flow and design variables converge. Performing one step involving unsteady computation of the flow field requires considerably less time than that required to solve the full steady-state equations via ordinal methods. Consequently, the iSA method can enable relatively fast computation. In practice, we solve the following time evolving equations: where M k is introduced to satisfy volume constraint x φ(x, t k+1 ) =V . The termination is considered by convergence of the design variable. If φ k+1 −φ k ≤ δ, then the algorithm terminates.

Numerical analysis of a liquid atomizer
The boundary conditions of computation is as follows. The inlet and outlet boundary conditions are defined on Γ in Fig. 3 Flow chart of the iSA method and Γ out in Fig. 2. Uniform velocityŪ is defined on inlet boundary, and pressure is fixed on outlet. The Reynolds number is 1.0×10 6 . The wall boundary condition is non-slip condition, where we used bounce-back scheme [43]. The position of the wall is calculated from the level-set function; if φ(x, t k ) < 0.5 then the point x is fluid, otherwise x is solid. The Smagorinsky model is used as the turbulence model.
The computation grid was determined considering the mesh resolution of the thinnest area of four channels. The thinnest area of the channel locates at the inlet to the swirl pool. The turbulence model we used is the LES mode, which needs a sufficient number of grids. Hence, we have divided the channel by 50 grids in the direction of channel diameter, and the number of grids was 800 × 800 × 800, corresponding to a total of 0.512 billion grids. In the numerical experiment in section 3, computation was carried out on 512 processors of XEON X5690. Consequently, 1 million grids are assigned to each processor.

Optimization of the atomizer
First, the initial shape is designed by an engineer. The part is designed to be fabricated using AM, as shown in Fig. 4(a), and is different from conventional shapes. The fuel comes from the bottom of the part to the liquid pool, and split into four channels. Then, the four channels are connected to the swirl chamber at the top of the part. The designer spent less than 1 day to design the initial shape. The designer set a wide cross-sectional area on the inlet of the flow channel, and put the inlet at the bottom of the pool aiming to gather fuel from the duct. The flow path is Fig. 4 Initial and optimal shapes. The smoothed shape is obtained by smoothing the optimal shape keeping the main flow path but smoothly decreasing the cross sectional area of the flow path then designed to smoothly connect to the swirl pool. The initial atomizer involved smoothly curved surfaces to ensure that the channels did not have any sharp edges, thereby leading to a low energy loss. The entrance of the channel was wide to ensure that the flow from the pool could enter the channels smoothly.
Then, the iSA method is applied, and the obtained solution is shown in Fig. 4(b). The optimization have terminated after 7500 steps. The computation time was about 66 h, which is acceptable from practical point of view. If the iSA method is not used, the computation time would be much larger. The connection area between the liquid pool and four channels was considerably reduced in the optimal shape compared to that in the initial shape. Moreover, the connection area was located around the middle of the duct, whereas it was located in the bottom area in the initial solution. The obtained optimal shape was represented by voxels, and thus, its surface was not smooth. Therefore, we manually set the median surface from the solution and generated a smoothed shape, as shown in Fig. 4(c). The smoothed shape retained the key feature of the optimal shape that the connection area was relatively small and located in the middle of the duct. The fuel comes from the bottom of the part to the middle of the liquid pool, and then split into four channels. The inlet areas of the four channels are relatively smaller than those of the initial shape. However, the area distribution along its stream line of the four channels is different from the obtained optimal shape; the area of the channel is gradually decreased along the stream line.

Verification through numerical analysis
The fluid pressure loss reduction was verified by performing a numerical analysis. Since conventional nozzle parts were in general designed using conventional computational fluid dynamics method, we used the same CFD method instead of the LBM approach, to analyze the initial and smoothed optimal shapes. Furthermore, a CFD analysis of the conventional shape (manufactured via machining tools) was performed for comparison. The conventional shape represents proprietary information and is not shown in this manuscript. Figure 5 shows the streamlines, and Table 1 shows the result of  the fluid energy loss. The energy loss of the initial and optimal shapes is considerably smaller than that of the conventional shape, because the channels of the initial and optimal shapes represent 3D-curved smooth paths owing to the manufacturing capability of AM, whereas the channels of the conventional shape involve right angles. Although the energy loss of the initial solution is significantly smaller than that of the conventional shape, the optimal shape indicated a considerably lower energy loss than the initial solution. The reason is understood by comparing the streamlines of the initial and optimal shapes. Specifically, the streamline of the junction between the pool and the channel is smoother in the optimal shape than that in the initial shape. The optimal flow path thus leads the liquid directly into the channels, whereas the flow of the initial shape meanders in the pool (Fig. 5(a)), leading to a higher energy loss.

Fabrication
We have fabricated the initial and optimal shapes using SLM 280, which is an AM machine by SLM solutions GmbH. The material was Ti-6Al-4V alloy. Stereolithography (STL) files were created based on the initial and optimal solutions. The print direction is from bottom to the top of the atomizer part aiming to conserve the axial symmetry after manufacturing. The channel consists of vertical surfaces and it is better to use support materials to obtain a smooth surface. However, support materials cannot be used because the supports in the channel cannot be removed after build. The half-cut test piece of the initial shape is shown in Fig. 6(a). It can be noted that the surface was rough, and it was required to be polished to obtain the smooth surface. The test pieces were subjected to abrasive flow machining (AFM) [27], and a smooth surface was obtained, as shown in Fig. 6(b) and (c). However, the shape can be noted to be modified after the application of AFM. Therefore, it is

Conclusion
A liquid atomizer was optimized by using the fluid topology optimization method, and fabricated by metal AM. Atomizers are generally manufactured using conventional machining tools, and their shapes are thus limited by the manufacturing capability. Because the shape of flow channels is vital to ensure the performance of fluid devices, it is desirable to combine topology optimization and AM approaches to design and manufacture fluid devices. The device fabricated in this study demonstrates that such a coupling is both feasible and applicable for industrial applications. Industrial applications usually needs large computation grid to compute flow fields, which leads to large computation time for optimization. The iSA method, which is a time-efficient topology optimization method, is employed to reduce computation time, and it output reasonable solution in a short time period. When designing industrial parts, it is believed that AM and topology optimization are essentially useful. However, the computation time is an inevitable issue in design process. How to modify the optimal shape is not studied. The present study shows that the iSA method is useful for optimizing the industrial parts. The total time needed to design the shape including initial shape design, optimization, and smoothing was about 3 days. When designing parts, especially aerospace parts, it took more than few weeks or even few months. Hence, 3 days is affordable in the industrial design.
However, parts fabricated using AM may exhibit surface roughness, thereby requiring further processing in the form of polishing. In the case of an optimized shape, polishing the fabricated part may change the shape of the part. Furthermore, the surface roughness is different among locations, materials, and characteristics of AM machines. Therefore, in future work, it is necessary to consider surface roughness in topology optimization.
Author contribution All authors contributed to the study conception and design. Optimization computation was performed by Kazuo Yonekura. CAD data handling and CFD computation were performed by Hitoshi Hattori. Fabrication was performed by Takafumi Nishizu. The first draft of the manuscript was written by Kazuo Yonekura and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript.
Funding Open access funding provided by The University of Tokyo. The authors declare that no funds, grants, or other support were received during the preparation of this manuscript.

Material availability
The data and material are not available.
Code availability The codes are not available.

Declarations
Consent for publication Informed consent for publication were obtained from all participants.

Competing interests
The authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons. org/licenses/by/4.0/.