1.1.2. Literature Search

In order to collect relevant literature for this review, a literature search was performed using Google Scholar based on the keyword combinations of "topology optimization" with the following:


In addition to the above, reverse tracking was used of citations of the seminal papers in the area, as well as relevant references in the papers from the search. Only journal publications have been included, except when important contributions have been made in available conference proceedings.

#### 1.1.3. Optimisation Methodology

A broad and open definition of *topology optimisation* is used herein. The presented methodologies must be capable of handling topological changes in three dimensions. That is, the methods should in a three-dimensional version be capable of handling large design changes and topological changes by *creating*, *removing* and *merging* holes. This can be difficult for some representations, especially in two dimensions, where auxiliary information such as topological derivatives is necessary for the creation of holes/structure. Pure shape optimisation, where only small modifications of the fluid–solid or fluid–void interface is possible, is not included in this review.

However, there might not be much need for changes in topology for most two-dimensional fluid flow problems. Due to the nature of fluid flow and the obvious objective of minimising power dissipation (or pressure drop), there is a desire to minimise the number of flow channels, i.e., only a single flow path is needed and the interface shape is modified. The need for topology optimisation does arise when other objectives, such as flow uniformity or diodicity, are considered and when the fluid flow is coupled to additional physics, such as e.g., heat transport.

The design representation used for topology optimisation of fluid-based problems is in general similar to those applied within the area of solid mechanics [2,6]. Figure 1 shows the three general options for representing the design. The first representation is an explicit boundary representation based on a body-fitted mesh adopting to the nominal geometry shown in red. If the design is changed, boundary nodes must be moved and the mesh must be updated or the domain must be entirely re-meshed if large changes are applied. The second representation is that of the density-based methods, which also includes level set methods where a smooth Heaviside projection is applied together with interpolation of material properties (so-called Ersatz material methods). The flow is penalised in the solid (black) domain, typically by modelling it as a porous material with very low permeability. The third representation is that of surface-capturing level set methods, where surface-capturing discretisation methods, e.g., the extended finite element method (X-FEM), where the cut elements are integrated using a special scheme and the interface boundary conditions are imposed, e.g., using stabilised Lagrange multipliers or a stabilised Nitsche's method.

**Figure 1.** Fluid nozzle illustrating the basic differences among design representations in topology optimisation: (**a**) explicit boundary representation (body fitted mesh); (**b**) density/ersatz material based representation; (**c**) level set based X-FEM/cutFEM representation.

The explicit boundary methods represent the physics well, but moving nodes and adaption of the mesh are non-smooth operations and this might pose difficulties in advancing the design for the optimiser. Furthermore, regularisation of the interface is necessary and, in case of full re-meshing, it might be difficult to assure high quality elements, while limiting the computational time.

The density-based methods are strong regarding the ability to change topology and change the design dramatically, due to the design sensitivities being distributed over a large part of the domain. The cost of introducing the design is relatively low, as only an extra term needs to be integrated, with no special interface treatment being necessary. However, there are problems such as choosing proper interpolation schemes for material properties and, in the case of fluids, a large enough penalisation of the flow in solid regions. The velocity and pressure fields are present in the entire domain, both solid and fluid regions, which may cause spurious flows and leaking pressure fields, if not penalised sufficiently.

For the surface-capturing methods, the well defined and crisp interface makes it easy to introduce interface couplings between different physics, e.g., for fluid–structure interaction problems. As for any level set method, due to the nature of the method, the design sensitivities are located only at the interface. This means that design changes can only propagate from the interface and no new holes appear automatically. This is often relieved by using an initial design with many holes or by introduction of a hole nucleation scheme, e.g., using topological derivatives.

The above methods will not be described in detail, but the readers are referred to the descriptions in the individual papers of the review and the general overview in the review papers by Sigmund and Maute [2] and van Dijk et al. [6].

#### *1.2. Layout of Paper*

The included papers are divided into different subsets based on the number and complexity of the physics involved. The layout of this paper is accordingly divided into sections.

The literature review is presented in Section 2. Section 2.1 covers pure fluid flow problems divided into steady laminar flow (Section 2.1.1), unsteady flow (Section 2.1.2), turbulent flow (Section 2.1.3) and non-Newtonian fluids (Section 2.1.4). Section 2.2 considers species transport problems. Section 2.3 deals with conjugate heat transfer problems, where the thermal field is modelled in both solid and fluid domains, divided according to the type of cooling into forced convection (Section 2.3.1) and natural convection (Section 2.3.2). Section 2.5 considers both material microstructures (Section 2.5.1), where effective material parameters are optimised, and porous media (Section 2.5.2), where homogenised properties are used to optimise a macroscale material distribution. Section 2.4 covers fluid–structure interaction (FSI) problems with the fluid flow loading a mechanical structure.

After the literature review, Section 3 performs a quantitative analysis of the included papers and Section 4 presents recommendations for future focus areas for research within topology optimisation of fluid-based problems. Finally, Section 5 briefly concludes the review paper.

#### **2. Literature Review**

In the following, the papers are grouped based on their most advanced example, if the work covers both simple and extended applications. Furthermore, the papers are presented in chronological order based on the date that the papers were available online, as this gives a better representation of the order than official date of the final issue, since that may well trail the online publication significantly and differently from paper to paper.

#### *2.1. Fluid Flow*

This section covers the majority of the papers included in the review paper, namely pure fluid flow problems. The section is divided into subsections covering steady laminar flow, unsteady flow, turbulent flow, microstructure and porous media.

#### 2.1.1. Steady Laminar Flow

Borrvall and Petersson [7] published the seminal work on fluid topology optimisation in 2003. They presented an in-depth mathematical basis for topology optimisation of Stokes flow. The design parametrisation is based on lubrication theory, leveraging the frictional resistance between parallel plates. Solid domains are approximated by areas with vanishing channel height. By designing the spatially-varying channel height, it is possible to achieve fluid topologies that dictate where flow channels minimising the dissipated energy are placed. This parametrisation was extended to

Navier–Stokes flow by Gersborg-Hansen et al. [8] in 2005. Both sets of authors note the similarity of the obtained equations with that of a Brinkman-type model of Darcy's law for flow through a porous medium. Here, solid domains are approximated by areas with a very low permeability. When treating two-dimensional problems of a finite depth, it makes sense to use the lubrication theory approach, since this ensures the out-of-plane viscous resistance due to finite channel width being taken into account in the fluid parts of the domain. However, for three-dimensional problems, the lubrication theory approach loses its physical meaning, whereas the porous media approach carries over without any issues.

Evgrafov [9] investigated the limits of porous materials in the topology optimisation of Stokes flow using mathematical analysis, complementing the analysis presented originally by Borrvall and Petersson [7]. Olesen et al. [10] presented a high-level programming-language implementation of topology optimisation for steady-state Navier–Stokes flow using the fictitious porous media approach. Guest and Prévost [11] took a different approach to the previous work and modelled the solid region as areas with Darcy flow of low permeability surrounded by areas of Stokes flow using an interpolated Darcy–Stokes finite element. Evgrafov [12] investigated the theoretical foundation and practical stability of the penalised Navier–Stokes equations going to the limit of infinite impermeability, showing that the problem is ill-posed for increasing impermeability of solid regions and that slight compressibility and filters can ensure solutions. Wiker et al. [13] treated problems with separate regions of Stokes and Darcy flow, similar to Guest and Prévost [11], but with a finite impermeability in order to simulate actual flow in porous media for a mass flow distribution problem. Pingen et al. [14] used the Lattice Boltzmann method (LBM) as an approximation of Navier–Stokes flow. Their work included the first three-dimensional result on a very coarse discretisation. Aage et al. [15] were the first to treat truly three-dimensional problems using shared-memory parallelisation, allowing them to optimise large scale Stokes flow problems. Bruns [16] highlighted the similarity of the previous approaches to that of a penalty formulation of imposing no-flow boundary conditions. Specifically, a volumetric penalty term is used to impose no-flow inside solid regions. This interpretation has become widespread through the years, as the physical relevance of the design parametrisation has lost interest and the strictly numerical view has gained popularity. Duan et al. [17–19] presented the first application of a variational level set method to fluid topology optimisation, producing similar designs to those previously obtained using density-based methods. Evgrafov et al. [20] performed a theoretical investigation into the use of kinetic theory to approximate Navier–Stokes for fluid topology optimisation. Othmer [21] derived a continuous adjoint formulation for both topology and shape optimisation of Navier–Stokes flow for implementation into the finite volume solver OpenFOAM. Although results for turbulent flow are presented, the approach is herein not considered "fully turbulent" due to the assumption of "frozen turbulence" (not taking the turbulence variables into account in the sensitivity analysis). Pingen et al. [22] discuss efficient methods for computation of the design sensitivities for LBM based methods and apply topology optimisation to a Tesla-type valve design problem. Zhou and Li [23] presented a variational level set method for Navier–Stokes flow excluding elements outside of the fluid region in order to increase accuracy of the no-slip boundary condition, but increasing book-keeping through updating the fluid mesh every design iteration. Pingen et al. [24] formulated a parametric level set approach using LBM to model the fluid flow. In contrast to the previous variational level set methods, gradient-based mathematical programming is used to update the level set function rather than using the traditional Hamilton–Jacobi equation. Similarly to Zhou and Li [23], Challis and Guest [25] proposed a variational level set method where only discrete fluid areas occur. This removes the need for interpolation schemes, and, by discarding the degrees-of-freedom in the solid, they are able to solve large three-dimensional problems at a reduced cost. Kreissl et al. [26] presented a generalised shape optimisation approach using an explicit level set formulation, where the nodal values of the level set are explicitly varied using a mathematical programming approach in contrast to variational [17–19,23,25] and parametric [24] level set formulations. Furthermore, they use a geometric boundary representation that enforces no-slip

conditions along the fluid–solid interface using the LBM. Liu et al. [27] optimised fluid distributors with multiple flow rate equality constraints. Okkels et al. [28] elaborated on the lubrication theory approach to design the height profile of the inlet of a bio-reactor minimising the pressure drop while maintaining a uniform flow through the reactor. While most of the previous works treated energy dissipation or pressure drop, Kondoh et al. [29] presented a density-based topology optimisation formulation for drag minimisation and lift maximisation of bodies using body force integration. Building on their previous work, Kreissl and Maute [30] developed an explicit level set method accurately capturing the boundary using the extended finite element method (X-FEM). Deng et al. [31] presented a density-based method for both steady and unsteady flows driven by gravitational, centrifugal and Coriolis body forces. Deng et al. [32] also proposed an implicit variational level set method for steady Navier–Stokes flow with body forces.

The work by Aage and Lazarov [33] is not focused on fluids, but rather presents a parallel framework for large-scale three-dimensional topology optimisation. However, they use a Stokes flow driven manifold distributor as an application using 3.35 million DOFs for the fluid problem. Evgrafov [34] presented a state space Newton's method for minimum dissipated energy of Stokes flow, which outperforms the traditional first order nested approaches significantly for specific problems. Romero and Silva [35] applied a density-based approach to the design of rotating laminar rotors by including the centrifugal and Coriolis forces. Yaji et al. [36] formulated a level set based topology optimisation of steady flow using the LBM and a continuous adjoint sensitivity analysis based on the Boltzmann equation. Liu et al. [37] presented a density-based method based on the LBM using a discrete adjoint sensitivity analysis posed in the moment space. Yonekura and Kanno [38] suggested to use transient information for steady state optimisation using the LBM by modifying the design during a single transient solve. Liu et al. [39] proposed a re-initialisation method for level set based optimisation using a regular mesh for the level set equation and an adaptive mesh for the fluid analysis. Duan et al. [40] presented an adaptive mesh method for density-based optimisation using an optimality criteria (OC) algorithm. Extending his previous work, Evgrafov [41] presented a Chebyshev method for topology optimisation of Stokes flow achieving locally cubic convergence for specific problems. Garcke and Hecht [42] performed an in-depth mathematical analysis of a phase field approach to topology optimisation of Stokes flow. Lin et al. [43] applied a density-based approach to the design of fixed-geometry fluid diodes, manufacturing the designs and verifying their performance using an experimental setup. Building on their previous formulation and analysis, Garcke et al. [44] provided numerical results for Stokes flow using their phase field approach and adaptive mesh refinement. Duan et al. [45] presented an Ersatz material based implicit level set method showing clear connections to density-based methods through the material distribution. Sá et al. [46] applied the concept of topological derivatives to fluid flow channel design.

Yoshimura et al. [47] proposed a gradient-free approach using a genetic algorithm to update a very coarse design using a Kriging surrogate model coupled to an immersed method known as the Building-Cube Method. Pereira et al. [48] applied a density-based approach to Stokes flow using polygonal elements and supplying a freely available code for fluid topology optimisation. Duan et al. [49] presented an adaptive mesh method to an Ersatz-material level set based approach for Navier–Stokes flow. Kubo et al. [50] used a variational Ersatz-material level set method to design manifolds for flow uniformity in microchannel reactors. Jang and Lee [51] maximise the acoustic transmission loss in a chamber muffler with a constraint on the reverse-flow power dissipation modelled using the Navier–Stokes equations. Koch et al. [52] proposed a method for automatic conversion from two-dimensional Ersatz-based level set topology optimisation to NURBS-based shape optimisation. Sato et al. [53] applied density-based topology optimisation to the design of no-moving parts fluid valves using a Pareto front exploration method. Sá et al. [54] further extended the work in [46] to rotating domains. Yonekura and Kanno [55] extended their previous approach for updating the design during the computation of an unsteady flow field to a level set based method. Dai et al. [56] presented a piecewise constant Ersatz-material level set method, which essentially reduces it to a

density-based method. Shen et al. [57] formulated a three-phase interpolation model in Darcy–Stokes flow for optimising fluid devices with Stokes flow, Darcy flow and solid domains. Deng et al. [58] used an Ersatz-material level set method to optimise two-phase flow, using with a phase field approach for the fluid–fluid interface. Garcke et al. [59] extended their phase field approach to Navier–Stokes flow, presenting both in-depth mathematical analysis and optimisation results for minimum drag and maximum lift problems. Alonso et al. [60] used a density-based approach to design laminar rotating swirl flow devices using a rotating axisymmetric model.

Jensen [61] proposed to use anisotropic mesh adaptation for density-based optimisation of Stokes flow, demonstrating that this allows for efficiently resolving the physical length scale related to very high Brinkman penalisation terms. Sá et al. [62] applied density-based optimisation to the design of a small scale rotary pump considering both energy dissipation and vorticity with experimental verification of the design performance. Zhou et al. [63] presented an integrated shape morphing and topology optimisation approach based on a mesh handling methodology called a deformable simplicial complex (DSC). Shin et al. [64] used a 2D axisymmetric model at moderate Reynolds numbers to optimise a vortex-type fluid diode that is postprocessed and verified by simulation under actual flow conditions and turbulence. Yonekura and Kanno [65] proposed a heuristic approximation of the Hessian matrix for fast density-based optimisation using the LBM. Behrou et al. [66] presented a methodology for adaptive explicit no-slip boundary conditions with a density-based method for laminar flow problems with mass flow constraints, adaptively removing elements in the solid regions. Alonso et al. [67] applied density-based optimisation to the design of Tesla-type centrifugal pumps without blades. Lim et al. [68] applied a density-based method to the design of vortex-type passive fluidic diode valves for nuclear applications. Sato et al. [69] presented a topology optimisation method for rarefied gas flow problems covering both gas and solid domains using the LBM. Gaymann et al. [70] applied a density-based method to the design of fluidic diode valves in both two- and three-dimensions at medium-to-high Reynolds numbers. Gaymann and Montomoli [71] applied deep neural networks and Monte Carlo Tree search to an absurdly coarse design grid.

#### 2.1.2. Unsteady Flow

All of the previous works consider steady-state flow problems. However, not all fluid flows develop a steady state and there are situations where the unsteady motion can be exploited by the design. This could either be due to transients in an onset flow or due to vortex shedding from an obstacle. It should be noted that some discretisation methods, e.g., lattice Boltzmann methods, are of a transient nature; however, the papers are categorised based on the use of a time-dependent objective functional and sensitivity analysis. Thus, if only the final state solution is used to update the design, it is considered steady state and not included in this section.

Kreissl et al. [72] presented the first work treating density-based topology optimisation for unsteady flow, using a discrete transient adjoint formulation and applied to the design of a diffuser and an oscillating flow manifold. They highlighted some difficulties encountered using a standard density approach combined with a stabilised finite element formulation. A few months later, Deng et al. [73] also presented density-based topology optimisation for unsteady flow, but using a continuous transient adjoint formulation. They optimised a wide arrangement of problems including oscillating manifold flows. Deng et al. [31] extended their work to both steady and unsteady flows driven by gravitational, centrifugal and Coriolis body forces. Abdelwahed and Hassine [74] introduced analysis and application of the topological gradient method for non-stationary fluid flows. Nørgaard et al. [75] presented topology optimisation of what can be considered the first truly unsteady flow problems, considering oscillating flow over multiple periods for both obstacle reconstruction and design of an oscillating pump using the LBM. Villanueva and Maute [76] formulated a CutFEM discretised explicit level set method to the design of two- and three-dimensional flow for both steady-state and fully transient problems. Although the use of a surface-capturing scheme, they show that there is still a penaltyand mesh-dependent mass loss through the interface. Chen et al. [77] presented a local-in-time

approximate discrete adjoint sensitivity analysis for unsteady flow using the LBM. They approximate the true time-dependent adjoint equations by splitting the time series into subintervals, rather than the full time series forward and then the full time series backwards. This significantly reduces the storage requirement but introduces an approximation error. Nørgaard et al. [78] discussed applications of automatic differentiation for topology optimisation, demonstrating their approach on an unsteady oscillating pressure pump. This was obtained using LBM and the design seems to rely on the slightly compressible behaviour of the fluid model. Sasaki et al. [79] presented fluid topology optimisation using a particle method for the first time. Specifically, they use the transient moving particle semi-implicit (MPS) method allowing them to treat free surface fluid flows without explicit surface tracking.
