Level set-based topology optimization for two dimensional turbulent flow using an immersed boundary method

Abstract

This paper presents a topology optimization method for two dimensional turbulent flow based on the Reynolds-averaged Navier-Stokes (RANS) equations using a level set boundary expression and the immersed boundary method (IBM). In this study, two-equation turbulence models, the k-ϵ and the k-ω, are considered. In our proposed method, the level set method is used for capturing the exact fluid-solid interfaces. Additionally, the no-slip boundary condition along the fluid-solid interfaces is imposed explicitly using the IBM during the optimization process. Based on the information of the exact boundary position, the interpolated velocity and pressure values within the viscous sublayer are estimated using the standard wall function. From the above formulations, we construct a topology optimization method for the total pressure drop minimization problems considering two dimensional turbulent flows, under the frozen turbulence assumption. We provide numerical examples to confirm the validity and utility of the proposed method.

Introduction

Currently, complex fluid flow in engineering applications can be analyzed using Computational Fluid Dynamics (CFD) through the use of commercial or open-source software. Additionally, due to steady progress in CPU performance, optimization approaches have become more attractive and realistic in engineering design applications of fluid flow.

With respect to optimization approaches, especially based on the adjoint method, after the pioneering work of Pironneau [1], wide-ranging research for fluid dynamics problems has been conducted [2], [3], [4], [5], [6]. A shape optimization approach is most often applied in aerodynamic optimization problems, such as the optimal shape design of NACA airfoils [7], turbine blades [8], and aircraft wings [6], since small modifications in shape have sufficient impact on the improvement of a performance, e.g., when drag minimization is an objective. In a shape optimization approach, the specific boundaries of the target object are changed to minimize or maximize a specific objective function based on the shape sensitivity, defined as the change in the objective functional in response to small perturbations of the structural boundaries.

In addition to research in shape optimization approaches, the topology optimization approach [9] is an attractive optimization method that offers certain advantages compared with sizing and shape optimization methods, due to its superior flexibility with respect to configuration changes, and allowance of the creation and disappearance of holes during the optimization, which can lead to innovative optimal configurations. The basic concepts of topology optimization are the extension of a design domain to a fixed design domain and the replacement of a structural optimization problem with a material distribution problem, using a characteristic function. The characteristic function used in the density approach, a widely-used topology optimization approach also termed the solid isotropic material with penalization (SIMP) method [10], represents a fictitious isotropic material whose elasticity tensor is assumed to be a function of penalized material density. Following this pioneering work, topology optimization methods have been widely applied to a variety of structural optimization problems such as stiffness maximization problems [11], [12], vibration problems [13], [14], optimum design problems for compliant mechanisms [15], [16], and thermal problems [17], [18].

Concerning topology optimization methods applied to fluid problems, Borrvall and Petersson first proposed a method for minimum power dissipation under Stokes flow [19], and this method was later extended to moderate Reynolds number (<1000) laminar flow regimes [20], [21], [22], [23], [24], [25], [26], [27], [28]. In these conventional approaches, employing the so-called Brinkman penalization method, the fixed design domain is assumed to be a porous medium by introducing the Darcy force term as an external source term in the Navier-Stokes equation, where the local porosities are considered as design variables.

From an engineering point of view, fluid behavior in almost all fluid applications can be assumed to be turbulent flow, except in microscale fluid machinery such as microelectromechanical systems (MEMS). However, there has been relatively little research conducted on topology optimization for turbulent flow regimes, compared with the prevalence of shape optimization approaches.

When considering turbulent flow using CFD, several approaches are available for solving the equation for the eddy viscosity using so-called turbulence models, among them, the Reynolds-averaged Navier-Stokes (RANS) simulation. Particularly in this study, we construct the topology optimization method for a fluid problem with the RANS simulation introducing the standard wall function, in order to reduce the time required for solving the governing equations.

The Spalart-Allmaras (S-A) turbulence model applied in previous studies concerning topology optimization problems [29], [30], [31] is an example of a one-equation RANS turbulence model typically used in external flow analyses, such as for aerodynamic problems. Papoutsis-Kiachagias et al. [29] addressed a topology optimization problem under incompressible laminar and turbulent ducted flows, especially for manifold designs, and also extended their proposed method to a heat transfer problem. In their approach, based on the Brinkman penalization method, new porosity-dependent terms were added to the main governing equations of the heat transfer and turbulence models. Kontoleontos et al. [30] extended a heat transfer problem and imposed constraints on the outlet flow direction, rates and mean outlet temperatures. In their design sensitivity analysis, the adjoint to the S-A turbulence model equation is taken into account for a continuous adjoint approach, and the local porosities as design variables are updated depending on the design sensitivity, using the steepest descent method. Yoon [31] explored topology optimization for turbulent flow using the S-A model and revealed the importance of eddy viscosity effects upon the optimal designs of several ducted flows. In this approach, the wall equation represented by the Eikonal equation was considered to calculate the distance value from the closest wall. Dilgen et al. [32] proposed the topology optimization method for turbulent flow considering one- and two-equation turbulence models without any simplifying assumptions in the sensitivity analysis. They compared the exact difference sensitivity with the design sensitivity derived from the sensitivity analysis under the frozen turbulence assumption, and revealed that differences between them are not negligible. All of these previous topology optimization approaches for turbulent flows used the Brinkman penalization method, so the no-slip boundary condition at the fluid-solid interface cannot be imposed explicitly during the optimization.

When considering the physical and numerical characteristics in turbulent flows, the imposition of the wall function requires clear boundaries in the expression of the fluid-solid interface in the topology optimization. However, in conventional topology optimization approaches using the Brinkman penalization method in fluid problems, the fluid-solid interface lacks clear boundaries because the interface is expressed as a porous medium. Also, when the Brinkman penalization method is used, the no-slip boundary condition cannot be explicitly applied to the fluid-solid interface because the porous medium by nature has an indistinct interface. Consequently, the fluid velocity and pressure distribution near the fluid-solid interface of the obtained optimized design is passively determined depending on the initial settings of an artificial numerical coefficient of the Darcy force term, where larger inverse permeability coefficient values can cause precipitous changes in flow velocity and pressure near the fluid-solid interface. This passive determination of the velocity and pressure distribution may cause unrealistic profiles of the velocity and pressure. These inexact phenomena can not be avoided even though the mesh refinement is applied. The imprecise fluid behavior treatments near the wall may therefore lead to unrealistic optimal designs, especially when dealing with topology optimization problems under turbulent flow.

In contrast to the above approaches, Challis and Guest presented a topology optimization method for Stokes flow that includes an explicitly enforced no-slip boundary condition during the optimization process [33]. They introduced the topological derivative as the design sensitivity and solved the Hamilton-Jacobi equation to evolve the level set function during the topology optimization. In a topology optimization problem for Navier-Stokes flow, Deng et al. [34] showed that the topological derivative can be considered as a weighted sum in the Hamilton-Jacobi equation for the evolution of the level set function. In their method, a reasonable weighting value for the topological derivative was determined after numerical examination prior to the optimization. This problematic need to set the values of certain parameters by hand also appears in other studies [35], [36], [37].

To fundamentally overcome the problem of the imprecise expression of the fluid-solid interface, there are two important requirements: 1) explicit imposition of the no-slip boundary condition on the interface during the topology optimization process, and 2) a method that expresses the fluid-solid interface with clear boundaries. The immersed boundary method (IBM) is a popular technique used in CFD analyses to explicitly impose the no-slip boundary condition. The IBM proposed by Peskin [38], [39] enforces the no-slip boundary condition by adding a body force to the Navier-Stokes equation as a reaction force from an object, using a fixed Cartesian grid. In particular, the IBM in a discrete forcing approach proposed by Fadlun et al. [40] is widely used, and several extensions of this method have been proposed [41], [42], [43], [44], [45]. In this method, the body force is introduced after the governing equations are discretized, which allows the no-slip boundary condition to be imposed more directly than in a continuous forcing approach in which the body force is incorporated into the governing equations before discretization. Several other IBM approaches are clearly reviewed and categorized by Mittal et al. [46].

In pioneering research in topology optimization for turbulent flow using the IBM, Sarstedt et al. [47] proposed a topology optimization method based on local optimality criteria using k-ϵ and k-ω SST turbulence models. In their optimization approach, local optimality criteria are used to deal with the minimization of the pressure loss under a specific volume constraint. They also provided the validity of the IBM in a discrete forcing approach with direct imposition of boundary conditions for the representation of the fluid-solid interface. From their result, corresponding to arbitrary velocity profiles in the fixed design domain, the IBM was in better agreement with the body-fitted mesh solution compared to the solution without using the IBM. Agreement was ensured even in $Re=2500$ turbulent flow, whereas differences increase with higher Re values unless the IBM is used.

For clear expression of the fluid-solid interface, the second critical requirement for precise evaluation of the fluid behavior near the wall under turbulent flow, a level set method is attractive. Following the basic methodology for tracking fronts and free boundaries proposed by Osher and Sethian [48], a level set method proposed for structural optimization [49], [12] enjoyed wide use for structural topology optimization problems, along with the SIMP method. Comprehensive literature concerning structural topology optimization based on level set methods is reviewed by Dijk et al. [50]. In this method, the design domain is fundamentally free of grayscales because structural boundaries are represented as the iso-surface of a scalar function, the level set function (LSF). Corresponding to a level set-based topology optimization for fluid problem, there are several recent studies [51], [52], [53].

In a previous study of a level set-based topology optimization method using the IBM, Kreissl et al. [54] proposed a topology optimization method for a laminar flow problem using an XFEM formulation of the incompressible Navier-Stokes equations, without employing the Brinkman approach. In their method, the no-slip boundary condition along the fluid-solid interface is enforced with a stabilized Lagrange multiplier method. Accordingly, the use of the XFEM inhibits the unrealistic flow penetration through thin walls that can occur when applying Brinkman penalization.

In this study, we propose a level set-based topology optimization method for ducted flows considering turbulent flow without using the Brinkman penalization approach to overcome the inexact distribution of the state variable of fluids at the fluid-solid interface. The turbulent flow is modeled using standard two-equation RANS models, i.e., k-ϵ and k-ω models, and these are discretized using the finite volume method (FVM).

In the proposed topology optimization process, we precisely and explicitly impose a no-slip boundary condition along the fluid-solid interface in the fixed design domain, using the IBM, a much different approach than the previous Brinkman penalization method in which the fluid-solid interface is expressed as a porosity. Corresponding to the governing equations, the Navier-Stokes equations are solved only in the fluid domain in the course of the our proposed topology optimization, while the solid domain is regarded as the porous medium and the Navier-Stokes equations are defined both in the fluid and solid domains in the Brinkman penalization approach. We implement the no-slip boundary condition explicitly on the fluid-solid interface of newly created holes during the topology optimization process and introduce the IBM of the discrete forcing approach, with direct imposition of boundary conditions [46], [55], [56], [57]. Furthermore, the LSF is used to obtain clear expressions of the fluid-solid interface. The calculation of distances from the wall required for wall function calculations is thereby facilitated.

In this study, we adopt the optimization algorithm in the level set-based topology optimization method proposed by Yamada et al. [58]. This method uses a topological derivative [59] as the design sensitivity and a reaction-diffusion equation as a time evolution equation of the LSF. A topological derivative is defined as a measure of the influence of creating an infinitesimally small hole in the fixed design domain, subject to appropriate boundary conditions. The advantages of this method are that not only are topological changes allowed during the optimization process, but a reinitialization treatment [60], typically required but time-consuming in conventional approaches to ensure accuracy when solving the Hamilton-Jacobi equation [61], [62], [63], is unnecessary.

The rest of this paper is as follows. Section 2 discusses the governing equations for turbulent flows, using the k-ϵ and k-ω turbulence models, and presents the wall function for the turbulent flow simulation. In Section 3, we formulate the topology optimization problem to minimize the total pressure drop in a steady-state incompressible viscous turbulent flow field, using a level set boundary expression. Based on this formulation, the topological derivative as the design sensitivity in this study is derived using the adjoint method, under the frozen turbulence assumption. In Section 4, the ghost-cell-based IBM of the discrete forcing approach, with direct imposition of boundary conditions, is presented and compared with the Brinkman penalization method. In Section 5, based on the optimization problem formulated in Section 3, we construct an optimization algorithm for ducted flows. In Section 6, we present numerical examples to verify the utility of the proposed optimization method, and provide a conclusion in Section 7.