Topology Optimization

Introduction

Topology optimization for fluid flow in AcuSolve is based on the solution of the flow equations combined with the solution of associated linearized adjoint equations. The objective is to minimize the mechanical losses of the fluid flow system in combination with constraints limiting the used material. To represent the solid geometry a field variable is used to indicate the presence of material which adds a resistance modeled by a porous media force with a positive Darcy coefficient. Refer to Sandboge et al. (2021) for details.

This document shows the equations used to solve the topology optimization problem for flow.

Optimization

Minimize mechanical losses
The considered optimization problems take the form:
Minimize
M ( u , γ ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcaaaa@3B73@
Subject to
D ( u , γ ) 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiraiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGHLjYScaaIWaaaaa@3DEA@
G k ( u , γ ) = 0 , k = 1 , ... , N MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4ramaaBa aaleaacaWGRbaabeaakiaacIcacaWG1bGaaiilaiabeo7aNjaacMca cqGH9aqpcaaIWaGaaiilaiaaykW7caaMc8UaaGPaVlaaykW7caaMc8 UaaGPaVlaaykW7caaMc8UaaGPaVlaaykW7caaMc8UaaGPaVlaaykW7 caaMc8UaaGPaVlaaykW7caaMc8UaaGPaVlaaykW7caaMc8UaaGPaVl aaykW7caaMc8UaaCzcaiaaxMaacaWGRbGaeyypa0JaaGymaiaacYca caGGUaGaaiOlaiaac6cacaGGSaGaamOtaaaa@6ABE@
R ( u , γ ) = 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGH9aqpcaaIWaaaaa@3D38@

where u=( u 1 , u 2 , u 3 ,p,...) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaaeaaaaaaaaa8 qacaWG1bGaeyypa0ZdaiaacIcapeGaamyDamaaBaaaleaacaaIXaaa beaakiaacYcacaWG1bWaaSbaaSqaaiaaikdaaeqaaOGaaiilaiaadw hadaWgaaWcbaGaaG4maaqabaGccaGGSaGaamiCaiaacYcacaGGUaGa aiOlaiaac6capaGaaiykaaaa@452D@ is the solution vector of flow variables u 1 MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamyDamaaBa aaleaacaaIXaaabeaaaaa@37D9@ for velocity vector component; p MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiCaaaa@36ED@ for pressure, and γ MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeq4SdCgaaa@379F@ is the design topology field which determines the presence of solid material. The function M(u,γ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcaaaa@3B73@ is the objective, D(u,γ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiraiaacI cacaWG1bGaaiilaiabeo7aNjaacMcaaaa@3B6A@ and G k (u,γ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4ramaaBa aaleaacaWGRbaabeaakiaacIcacaWG1bGaaiilaiabeo7aNjaacMca aaa@3C93@ are constraint functions. The flow equations are represented by R(u,γ)=0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGH9aqpcaaIWaaaaa@3D38@ acting as an additional constraint.

Here, M ( u , γ ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcaaaa@3B73@ represents the mechanical losses as a single objective, integrating over the inlet and outlet surfaces.

M ( u , γ ) = S I p + 1 2 ρ u k 2 u i n i d S S O p + 1 2 ρ u k 2 u i n i d S MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGH9aqpcqGHsisldaWdrbqa amaabmaabaGaamiCaiabgUcaRmaalaaabaGaaGymaaqaaiaaikdaaa GaeqyWdiNaamyDamaaDaaaleaacaWGRbaabaGaaGOmaaaaaOGaayjk aiaawMcaaaWcbaGaam4uamaaBaaameaacaWGjbaabeaaaSqab0Gaey 4kIipakiaaykW7caWG1bWaaSbaaSqaaiaadMgaaeqaaOGaaGPaVlaa d6gadaWgaaWcbaGaamyAaaqabaGccaaMc8UaamizaiaadofacqGHsi sldaWdrbqaamaabmaabaGaamiCaiabgUcaRmaalaaabaGaaGymaaqa aiaaikdaaaGaeqyWdiNaamyDamaaDaaaleaacaWGRbaabaGaaGOmaa aaaOGaayjkaiaawMcaaaWcbaGaam4uamaaBaaameaacaWGpbaabeaa aSqab0Gaey4kIipakiaaykW7caWG1bWaaSbaaSqaaiaadMgaaeqaaO GaaGPaVlaad6gadaWgaaWcbaGaamyAaaqabaGccaaMc8Uaamizaiaa dofaaaa@6EF3@

where ρ MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdihaaa@37B8@ is the density, S MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4uaaaa@36D0@ is the surface.

The function D ( u , γ ) = D ( γ ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiraiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGH9aqpcaWGebGaaiikaiab eo7aNjaacMcaaaa@4039@ is here the function of design fraction occupying the solid volume and is used as an optional constraint to limit the amount of solid material.

D ( γ ) = Ω γ d Ω Ω d Ω MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiraiaacI cacqaHZoWzcaGGPaGaeyypa0ZaaSaaaeaadaWdrbqaaiabeo7aNjaa ykW7caWGKbGaeyyQdCfaleaacqGHPoWvaeqaniabgUIiYdaakeaada WdrbqaaiaadsgacqGHPoWvaSqaaiabgM6axbqab0Gaey4kIipaaaaa aa@4AB4@

where Ω MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeyyQdCfaaa@3787@ is the volume.

As an alternative, the constraint can also be given as an upper bound constraint D(u,γ)0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiraiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGHKjYOcaaIWaaaaa@3DD9@ .

Further, the functions G k ( u , γ ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4ramaaBa aaleaacaWGRbaabeaakiaacIcacaWG1bGaaiilaiabeo7aNjaacMca aaa@3C93@ are used to constrain mass flow rates at outlets 𝑘 for the 𝑁 outlets,

G k ( u , γ ) = 1 2 m l m ¯ l 2 = 1 2 S O k ρ u j n j d S m ¯ k 2 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4ramaaBa aaleaacaWGRbaabeaakiaacIcacaWG1bGaaiilaiabeo7aNjaacMca cqGH9aqpdaWcaaqaaiaaigdaaeaacaaIYaaaamaabmaabaGaamyBam aaBaaaleaacaWGSbaabeaakiabgkHiTiqad2gagaqeamaaBaaaleaa caWGSbaabeaaaOGaayjkaiaawMcaamaaCaaaleqabaGaaGOmaaaaki abg2da9maalaaabaGaaGymaaqaaiaaikdaaaWaaeWaaeaadaWdrbqa aiabeg8aYjaadwhadaWgaaWcbaGaamOAaaqabaGccaWGUbWaaSbaaS qaaiaadQgaaeqaaOGaaGPaVlaadsgacaWGtbaaleaacaWGtbWaaSba aWqaaiaad+eadaWgaaqaaiaadUgaaeqaaaqabaaaleqaniabgUIiYd GccqGHsislceWGTbGbaebadaWgaaWcbaGaam4AaaqabaaakiaawIca caGLPaaadaahaaWcbeqaaiaaikdaaaaaaa@5D7B@

where m k MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamyBamaaBa aaleaacaWGRbaabeaaaaa@3806@ is the mass flow rate at exit 𝑘 and m ¯ k MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGabmyBayaara WaaSbaaSqaaiaadUgaaeqaaaaa@381E@ is the target mass flow rate at the same exit.

The constraints are naturally optional. It is possible to have an optimization problem without any constraints, except for the flow constraint Equation 4 which is always active.

Minimizing the mechanical losses and fluid volume
If M ( u , γ ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcaaaa@3B73@ represents the mechanical losses in Equation 1, regions of zero velocity will not change the objective value. Thus, zero velocity "pockets" in the solid region can occur in the final solution. To avoid such unclean solid geometry, it is possible to define M ( u , γ ) MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcaaaa@3B73@ as a multi-objective adding the size of the fluid volume in the design domain to the mechanical losses. The strength of minimization of the fluid volume size is defined by a user-defined coefficient c s MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4yamaaBa aaleaacaWGZbaabeaaaaa@3803@ multiplied by a normalization coefficient ε n MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyTdu2aaS baaSqaaiaad6gaaeqaaaaa@38BD@ ,
Minimize
M ( u , γ ) + c s ε n D γ MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGHRaWkcaWGJbWaaSbaaSqa aiaadohaaeqaaOGaeqyTdu2aaSbaaSqaaiaad6gaaeqaaOGaamiram aabmaabaGaeq4SdCgacaGLOaGaayzkaaaaaa@4534@
Subject to
R(u,γ)=0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGH9aqpcaaIWaaaaa@3D38@

Flow Equations

Computational domain
The flow equations are solved in a computational domain which is allowed to vary during the optimization process by modifying the material properties. At the first iteration, the design domain may use material properties representing a fluid, such as air, see Figure 1, and at convergence some parts of the domain have built up regions with material properties representing a solid material, see Figures Figure 2Figure 3. The computational domain Ω MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeyyQdCfaaa@3787@ is fixed, but the material properties are functions of a design field 𝛾 as is explained in the continuous adjoint for topology optimization section.
Figure 1. Design domain with all fluid material properties. A starting condition with sub-optimal mechanical losses.


Figure 2. Topology solution where some parts of the design domain have material properties representing a solid material in order to minimize the mechanical losses of the system.


Figure 3. Topology solution minimizing the mechanical losses subject to unequal mass flow constraints at the outlets.


Incompressible Navier-Stokes equations

The state equation residuals of the incompressible equations are defined as, R u i = 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuamaaBa aaleaacaWG1bGaamyAaaqabaGccqGH9aqpcaaIWaaaaa@3AAC@ , R p = 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuamaaBa aaleaacaWGWbaabeaakiabg2da9iaaicdaaaa@39B9@ where

R u i = ρ u i t + ρ u j u i x j + p x j τ i j x j ρ b ^ i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuamaaBa aaleaacaWG1bGaamyAaaqabaGccqGH9aqpcqaHbpGCdaWcaaqaaiab gkGi2kaadwhadaWgaaWcbaGaamyAaaqabaaakeaacqGHciITcaWG0b aaaiabgUcaRiabeg8aYjaadwhadaWgaaWcbaGaamOAaaqabaGcdaWc aaqaaiabgkGi2kaadwhadaWgaaWcbaGaamyAaaqabaaakeaacqGHci ITcaWG4bWaaSbaaSqaaiaadQgaaeqaaaaakiabgUcaRmaalaaabaGa eyOaIyRaamiCaaqaaiabgkGi2kaadIhadaWgaaWcbaGaamOAaaqaba aaaOGaeyOeI0YaaSaaaeaacqGHciITcqaHepaDdaWgaaWcbaGaamyA aiaadQgaaeqaaaGcbaGaeyOaIyRaamiEamaaBaaaleaacaWGQbaabe aaaaGccqGHsislcqaHbpGCceWGIbGbaKaadaWgaaWcbaGaamyAaaqa baaaaa@62D8@
R p = u i x i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuamaaBa aaleaacaWGWbaabeaakiabg2da9maalaaabaGaeyOaIyRaamyDamaa BaaaleaacaWGPbaabeaaaOqaaiabgkGi2kaadIhadaWgaaWcbaGaam yAaaqabaaaaaaa@4010@

where p MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiCaaaa@36EC@ , u i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamyDamaaBa aaleaacaWGPbaabeaaaaa@380B@ denote the static pressure and velocity components, respectively and 𝜌 denotes density and b ^ i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGabmOyayaaja WaaSbaaSqaaiaadMgaaeqaaaaa@3808@ are the specific body force components. The stresses τ i j MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqiXdq3aaS baaSqaaiaadMgacaWGQbaabeaaaaa@39C5@ are given by

τ ij =μ u i x j + u j x i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqiXdq3aaS baaSqaaiaadMgacaWGQbaabeaakiabg2da9iabeY7aTnaabmaabaWa aSaaaeaacqGHciITcaWG1bWaaSbaaSqaaiaadMgaaeqaaaGcbaGaey OaIyRaamiEamaaBaaaleaacaWGQbaabeaaaaGccqGHRaWkdaWcaaqa aiabgkGi2kaadwhadaWgaaWcbaGaamOAaaqabaaakeaacqGHciITca WG4bWaaSbaaSqaaiaadMgaaeqaaaaaaOGaayjkaiaawMcaaaaa@4D2E@

where μ MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqiVd0gaaa@37AD@ is the viscosity.

Continuous Adjoint for Topology Optimization

The optimization problem (Equation 1) is solved using the adjoint methods. The continuous adjoint approach for topology optimization is presented in this section.

Momentum equations for topology optimization

The Darcy porosity model for the momentum flow equations is used; the residual is expressed as

R u i p = ρ u i t + ρ u j u i x j + p x j τ i j x j ρ b ^ i + f i = 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuamaaDa aaleaacaWG1bGaamyAaaqaaiaadchaaaGccqGH9aqpcqaHbpGCdaWc aaqaaiabgkGi2kaadwhadaWgaaWcbaGaamyAaaqabaaakeaacqGHci ITcaWG0baaaiabgUcaRiabeg8aYjaadwhadaWgaaWcbaGaamOAaaqa baGcdaWcaaqaaiabgkGi2kaadwhadaWgaaWcbaGaamyAaaqabaaake aacqGHciITcaWG4bWaaSbaaSqaaiaadQgaaeqaaaaakiabgUcaRmaa laaabaGaeyOaIyRaamiCaaqaaiabgkGi2kaadIhadaWgaaWcbaGaam OAaaqabaaaaOGaeyOeI0YaaSaaaeaacqGHciITcqaHepaDdaWgaaWc baGaamyAaiaadQgaaeqaaaGcbaGaeyOaIyRaamiEamaaBaaaleaaca WGQbaabeaaaaGccqGHsislcqaHbpGCceWGIbGbaKaadaWgaaWcbaGa amyAaaqabaGccqGHRaWkcaWGMbWaaSbaaSqaaiaadMgaaeqaaOGaey ypa0JaaGimaaaa@6889@

where f i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOzamaaBa aaleaacaWGPbaabeaaaaa@37FC@ are the components of the porous media force vector

f i = C Darcy γ μ u i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOzamaaBa aaleaacaWGPbaabeaakiabg2da9iaadoeadaWgaaWcbaGaaeiraiaa bggacaqGYbGaae4yaiaabMhaaeqaaOWaaeWaaeaacqaHZoWzaiaawI cacaGLPaaacaaMc8UaaGPaVlabeY7aTjaaykW7caWG1bWaaSbaaSqa aiaadMgaaeqaaaaa@4A27@

and C Darcy MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4qamaaBa aaleaacaqGebGaaeyyaiaabkhacaqGJbGaaeyEaaqabaaaaa@3B6D@ is the Darcy coefficient,

C Darcy γ = C Darcy f + γ p C Darcy s C Darcy f MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4qamaaBa aaleaacaqGebGaaeyyaiaabkhacaqGJbGaaeyEaaqabaGcdaqadaqa aiabeo7aNbGaayjkaiaawMcaaiabg2da9iaadoeadaqhaaWcbaGaae iraiaabggacaqGYbGaae4yaiaabMhaaeaacaWGMbaaaOGaey4kaSIa eq4SdC2aaWbaaSqabeaacaWGWbaaaOWaaeWaaeaacaWGdbWaa0baaS qaaiaabseacaqGHbGaaeOCaiaabogacaqG5baabaGaam4Caaaakiab gkHiTiaadoeadaqhaaWcbaGaaeiraiaabggacaqGYbGaae4yaiaabM haaeaacaWGMbaaaaGccaGLOaGaayzkaaaaaa@5929@

where the Solid Isotropic Material with Penalization (SIMP) defined by p MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiCaaaa@36EC@ , a constant whose value is selected p = 1 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiCaiabg2 da9iaaigdaaaa@38AD@ to p = 3 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiCaiabg2 da9iaaiodaaaa@38AF@ . Thus, the design field γ MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeq4SdCgaaa@379E@ determines the level of resistance in the design domain; zero resistance for fluid regions and a high resistance to minimize the flow through the solid regions.
Note: For γ = 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeq4SdCMaey ypa0JaaGimaaaa@395E@ you have C Darcy = 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4qamaaBa aaleaacaqGebGaaeyyaiaabkhacaqGJbGaaeyEaaqabaGccqGH9aqp caaIWaaaaa@3D37@ and Equation 13 becomes identical to Equation 10, and for γ = 1 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeq4SdCMaey ypa0JaaGymaaaa@395F@ you have C Darcy = C Darcy s MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4qamaaBa aaleaacaqGebGaaeyyaiaabkhacaqGJbGaaeyEaaqabaGccqGH9aqp caWGdbWaa0baaSqaaiaabseacaqGHbGaaeOCaiaabogacaqG5baaba Gaam4Caaaaaaa@42EC@ .

The value of C Darcy s MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4qamaaDa aaleaacaqGebGaaeyyaiaabkhacaqGJbGaaeyEaaqaaiaadohaaaaa aa@3C66@ is automatically computed by the optimization algorithm, aiming to minimize the flow through the solid region but at the same time limiting the size of the force f i MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOzamaaBa aaleaacaWGPbaabeaaaaa@37FC@ to avoid poor convergence.

Adjoint equations

The adjoint equations to the linearized incompressible Navier-Stokes equations are given by,

R u i a = ρ u j u i a x j + u j a x i + p a x j τ i j a x j + f i a + F Ω = 0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuamaaBa aaleaacaWG1bWaa0baaWqaaiaadMgaaeaacaWGHbaaaaWcbeaakiab g2da9iabgkHiTiabeg8aYjaadwhadaWgaaWcbaGaamOAaaqabaGcda qadaqaamaalaaabaGaeyOaIyRaamyDamaaDaaaleaacaWGPbaabaGa amyyaaaaaOqaaiabgkGi2kaadIhadaWgaaWcbaGaamOAaaqabaaaaO Gaey4kaSYaaSaaaeaacqGHciITcaWG1bWaa0baaSqaaiaadQgaaeaa caWGHbaaaaGcbaGaeyOaIyRaamiEamaaBaaaleaacaWGPbaabeaaaa aakiaawIcacaGLPaaacqGHRaWkdaWcaaqaaiabgkGi2kaadchadaah aaWcbeqaaiaadggaaaaakeaacqGHciITcaWG4bWaaSbaaSqaaiaadQ gaaeqaaaaakiabgkHiTmaalaaabaGaeyOaIyRaeqiXdq3aa0baaSqa aiaadMgacaWGQbaabaGaamyyaaaaaOqaaiabgkGi2kaadIhadaWgaa WcbaGaamOAaaqabaaaaOGaey4kaSIaamOzamaaDaaaleaacaWGPbaa baGaamyyaaaakiabgUcaRiaadAeadaWgaaWcbaGaeuyQdCfabeaaki abg2da9iaaicdaaaa@6DF3@
R p a = u i a x i =0 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuamaaBa aaleaacaWGWbWaaWbaaWqabeaacaWGHbaaaaWcbeaakiabg2da9maa laaabaGaeyOaIyRaamyDamaaDaaaleaacaWGPbaabaGaamyyaaaaaO qaaiabgkGi2kaadIhadaWgaaWcbaGaamyAaaqabaaaaOGaeyypa0Ja aGimaaaa@43E0@

where the adjoint variables u a MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamyDamaaCa aaleqabaGaamyyaaaaaaa@3804@ and p a MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiCamaaCa aaleqabaGaamyyaaaaaaa@37FF@ measures the sensitivities of varying the field γ MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeq4SdCgaaa@379E@ around the fluid solution u MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamyDaaaa@36F1@ and p MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamiCaaaa@36EC@ of the forward problem.

Optimization Algorithm

Derivative of responses

The variation of the mechanical losses objective function (Equation 5) with respect to the design field variables is given by the expression

F ˜ = Ω C ˜ Darcy γ μ k i u i u i a d Ω MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGabmOrayaaia Gaeyypa0Zaa8quaeaaceWGdbGbaGaadaWgaaWcbaGaaeiraiaabgga caqGYbGaae4yaiaabMhaaeqaaOWaaeWaaeaacqaHZoWzaiaawIcaca GLPaaaaSqaaiabgM6axbqab0Gaey4kIipakiaaykW7caaMc8UaaGPa VpaalaaabaGaeqiVd0gabaGaam4AamaaBaaaleaacaWGPbaabeaaaa GccaWG1bWaaSbaaSqaaiaadMgaaeqaaOGaamyDamaaDaaaleaacaWG PbaabaGaamyyaaaakiaaykW7caaMc8UaamizaiabgM6axbaa@57A6@

where u i a MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamyDamaaDa aaleaacaWGPbaabaGaamyyaaaaaaa@38F2@ is the adjoint velocity vector. The coefficient C ˜ Darcy γ MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGabm4qayaaia WaaSbaaSqaaiaabseacaqGHbGaaeOCaiaabogacaqG5baabeaakmaa bmaabaGaeq4SdCgacaGLOaGaayzkaaaaaa@3EB6@ is the variation of the darcy factor, which is constant in the design domain for the SIMP-method of order 1. In essence, the derivative of the design material is determined from the scalar product between the velocity and the adjoint velocity vectors.

A volume constraint based on the topology design fraction in Equation 6 gives the variation

D ˜ = D γ γ ˜ Ω d Ω = Ω γ ˜ d Ω Ω d Ω MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGabmirayaaia Gaeyypa0JabmirayaafaWaaeWaaeaacqaHZoWzaiaawIcacaGLPaaa daWcaaqaaiqbeo7aNzaaiaaabaWaa8quaeaacaWGKbGaeyyQdCfale aacqGHPoWvaeqaniabgUIiYdaaaOGaeyypa0ZaaSaaaeaadaWdrbqa aiqbeo7aNzaaiaGaaGPaVlaadsgacqGHPoWvaSqaaiabgM6axbqab0 Gaey4kIipaaOqaamaapefabaGaamizaiabgM6axbWcbaGaeyyQdCfa beqdcqGHRiI8aaaaaaa@54FF@

which gives a constant value in the design domain for the same volume variation γ ˜ MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGafq4SdCMbaG aaaaa@37AE@ .

The function G k u , γ MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4ramaaBa aaleaacaWGRbaabeaakmaabmaabaGaamyDaiaacYcacqaHZoWzaiaa wIcacaGLPaaaaaa@3CC4@ given by Equation 7, used to account for the constraint on mass flow at outlet k MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4Aaaaa@36E8@ , is specified as an equality constraint and actually treated as a double inequality constraint defined below as

δ k Tol m k m ¯ k MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqiTdq2aa0 baaSqaaiaadUgaaeaacaqGubGaae4BaiaabYgaaaGccqGHKjYOcaWG TbWaaSbaaSqaaiaadUgaaeqaaOGaeyOeI0IabmyBayaaraWaaSbaaS qaaiaadUgaaeqaaaaa@425C@

or

G k ( u , γ ) 1 2 δ k Tol 2 MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaam4ramaaBa aaleaacaWGRbaabeaakiaacIcacaWG1bGaaiilaiabeo7aNjaacMca cqGHKjYOdaWcaaqaaiaaigdaaeaacaaIYaaaamaabmaabaGaeqiTdq 2aa0baaSqaaiaadUgaaeaacaqGubGaae4BaiaabYgaaaaakiaawIca caGLPaaadaahaaWcbeqaaiaaikdaaaaaaa@47C5@

where δ k Tol MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqiTdq2aa0 baaSqaaiaadUgaaeaacaqGubGaae4BaiaabYgaaaaaaa@3B72@ is a user-defined tolerance value for the exit 𝑘. The variation of this function is equal to

G ˜ k = m k m ¯ k S O k ρ u ˜ j n j d S MathType@MTEF@5@5@+= feaahGart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGabm4rayaaia WaaSbaaSqaaiaadUgaaeqaaOGaeyypa0ZaaeWaaeaacaWGTbWaaSba aSqaaiaadUgaaeqaaOGaeyOeI0IabmyBayaaraWaaSbaaSqaaiaadU gaaeqaaaGccaGLOaGaayzkaaWaa8quaeaacqaHbpGCceWG1bGbaGaa daWgaaWcbaGaamOAaaqabaGccaWGUbWaaSbaaSqaaiaadQgaaeqaaO GaamizaiaadofaaSqaaiaadofadaWgaaadbaGaam4taiaadUgaaeqa aaWcbeqdcqGHRiI8aaaa@4CA1@

where 𝑆𝑂𝑘 is the surface for outlet 𝑘.

Given the values of the objectives and constraints, and the field derivatives with respect to 𝛾 for the objective, (Equation 18 for mechanical losses), and constraints, (constant for the design fraction constraint), the optimizer can update the design field 𝛾. The default optimizer is the methods of moving asymptotes (MMA). Refer to Svanberg (2007) for details of the MMA.

Metrics

The design variable 𝛾 is continuous, but a well-defined geometry should have the value either 𝛾 = 0 for a fluid or 𝛾 = 1 for a solid with small regions of the design domain where 0 < 𝛾 < 1. To measure the amount a well-defined solution of either fluid or solid, a quality index can be used,

I Q = 1 4 V γ 1 γ d V V MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamysamaaBa aaleaacaWGrbaabeaakiabg2da9iaaigdacqGHsislcaaI0aGaaGPa VlaaykW7daWcaaqaamaapefabaGaeq4SdC2aaeWaaeaacaaIXaGaey OeI0Iaeq4SdCgacaGLOaGaayzkaaaaleaacaWGwbaabeqdcqGHRiI8 aOGaaGPaVlaadsgacaWGwbaabaGaamOvaaaaaaa@4C3C@

where 0 ≤ 𝐼𝑄 ≤ 1, and 𝐼𝑄 = 1 for a well-defined geometry. Here, 𝑉 is the volume of the design domain.

Illustration of the optimization algorithm

Equation 14 provides information on how the design material should be altered in the design space at each design iteration. To determine if solid material should be added or removed, it is possible to view the velocity 𝑢 and adjoint velocity 𝑢𝑎 vectors in conjunction; if the angle is less than 90 degrees, material should be removed and if the angle is greater than 90 degrees material should be added.

Figure 4 shows (i) the initial velocity and adjoint velocity fields where the vectors have an angle greater than 90 degrees in some regions in front and rear of the obstacle; and (ii) the final velocity and adjoint velocity fields where all angles are less than 90 degrees in the entire fluid design domain.
Figure 4. Flow past an obstacle. The velocity fields (white) and adjoint velocity fields (red) at the initial and final iterations: top for the initial field, bottom for the final field.


A close-up of the region behind the obstacle is shown in Figure 5, where the velocity shows some circulation while the adjoint velocity is attached without circulation.
Figure 5. The velocity and adjoint velocity fields behind the obstacle.


Examples

Duct with sharp turn
As a simple example consider the duct domain (Figure 6), where the following optimization problem is solved,
Minimize
M ( u , γ ) MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcaaaa@3B74@
Subject to
R ( u , γ ) = 0 MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGH9aqpcaaIWaaaaa@3D39@
where 𝑀(𝑢,𝛾) is the function of Equation 5, representing the mechanical losses of the fluid flow.
Figure 6. Flow optimization in a duct.. Left: shows the empty design domain; right shows the topology solution of Equation 24 where the solid geometry is shown in grey.


To solve the multi objective problem where you also factor in the size of the occupied fluid volume you have,
Minimize
M ( u , γ ) + c s ε n D ( γ ) MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamytaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGHRaWkcaWGJbWaaSbaaSqa aiaadohaaeqaaOGaeqyTdu2aaSbaaSqaaiaad6gaaeqaaOGaamirai aacIcacqaHZoWzcaGGPaaaaa@4505@
Subject to
R ( u , γ ) = 0 MathType@MTEF@5@5@+= feaahGart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOuaiaacI cacaWG1bGaaiilaiabeo7aNjaacMcacqGH9aqpcaaIWaaaaa@3D39@
By using a larger coefficient 𝑐𝑠, the duct cross-section will tend to be smaller due to a higher weight minimizing the occupied fluid volume ().
Figure 7. Flow optimization in a duct.. Using Equation 26. Left; shows the topology solution using 𝑐𝑠 =1. Right; shows the topology solution using 𝑐𝑠 =25, here illustrated by the surface between the fluid and solid regions.


Design-optimization of an HVAC Ducting System
The optimization technology is illustrated in this section in the design of the HVAC system of the Altair CX1 car. The cabin geometry of the Altair car is illustrated in Figure 8.
Figure 8. Inlets to the cabin geometry and air supply locations from the blower.


In this figure the ten different inlets to the cabin are shown. The six of them are aiming at the front and rear upper part of the cabin and the remaining four are providing the air to the front and rear lower part of the cabin. The supply of the air from the blower has also been separated into two parts. The first supply has to be connected with the upper inlets and the second supply has to be connected to the lower inlets to the cabin.

The first part of the algorithmic procedure is based on topology optimization for minimum power dissipation or mechanical energy losses to account for minimum fuel consumption and minimum volume to account for minimum manufacturing cost constrained by specific mass flow rates at each exit. In this study all mass flow rates selected as constraints were equal to 0.001kg/s and the tolerance of the constraints were selected equal to 5%.

Given the available space from the cabin and generally the car components, the first step is to define the initial box within which the duct must be confined. Apart from the availability in space the only other restriction to apply for topology optimization is to define the inlet to the box as the air supply location from the blower(s) and the outlet from the box as the inlets to the cabin. The CAD geometry can be in general automatically created and can be extracted from the available space. In this case a set of simple boxes were created to define the initial domain (Figure 9).
Figure 9. Initial CAD domain based on the available space from the cabin.


The next step is to generate the computational mesh at the initial domain. This mesh can consist of either hexaedral or tetrehedral elements. For this study, a hexahedral mesh with 2.5 ×105 nodes is used, Figure 10, left, which has sufficient resolution for topology optimization.

The topology optimization of the two ducts was performed in a sequential manner. First, the duct connecting the one inlet with the six upper outlets (cabin inlets) was optimized. For that reason, a part of the domain was excluded from the solution at this step to allow the evolution of the design of the first duct in a way that would not block the evolution of the other duct. This exclusion was made possible by applying nodal boundary conditions on the design variable at the excluded part of the domain and not through separate meshing.

After, the first optimal duct geometry is computed, the space it requires is excluded from the initial box, using again nodal boundary conditions on the design variable. To avoid intersections and very close proximity between the two duct geometries a user selected area that surrounds the first duct is also excluded. The second duct is optimized, and the optimal geometry uses the same computational mesh. Thus, a single initial mesh should only be created and no remeshing is required while performing the optimization of the separate ducts. This procedure can be continued to a next round of optimization of the two ducts; the first duct can be optimized again, leaving out not the initial random space but the space covered by the second duct (plus an additional layer), then the second duct can be optimized again based on the new shape of the first duct and so forth. This approach is expected to result in designs with even better performance but was omitted in this study.

The final design of the two ducts using topology optimization is afterwards extracted and smoothed. The final shape is shown in Figure 10, right.
Figure 10. Computational mesh of the initial domain and optimal design of the HVAC ducts based on topology optimization.


A body fitted mesh with 6.5×105nodes is generated at this geometry in order to further optimize the shape using the shape optimization methodology with morph shapes. Shape optimization is used to further improve the value of the power dissipation objective function and ensure that the constraints on mass flow rates at the exits of the domain are as accurately as possible satisfied. It has been observed that the values of the mass flow rates at the exits are very sensitive to the geometrical design and mesh and, thus, the smoothing and remeshing of the optimized geometry affects them quite a lot. The body fitted mesh is also expected to give more accurate estimation of the quantities of interest and, thus shape optimization is expected to give the completely final design of the ducts. It should also be mentioned that the two ducts are handled as one during this step to avoid intersections and high proximity in the final design.

The morph shapes, required for the shape optimization described in the theory, are generated automatically using the flow field computed in the body fitted mesh. This way the two methodologies, topology and shape optimization described above are not disconnected but rather constitute the ingredients of a unified methodology. The flow velocity is the metric used to define optimally the morph shapes so that several types of morph shape can be formed as well as to ensure that the whole domain is appropriately space-filled and all the parts of the design are taken into consideration while performing the shape optimization step.

The automatically computed centers of the morph shapes are shown in Figure 11, left. In this case, the length scale of morph shape generator was equal to 0.2m which resulted in the generation of 22 control points in total. In Figure 11, right, the morph shape vectors are shown.
Figure 11. Bending and area morph shape centers and vectors.


Two different types of morph shapes are generated; the first one is depicted by the red vectors, and it allows for the bending of the duct during optimization. Two design variables per node are needed for this type of morph shape to cover all the design space. The second type of morph shapes is illustrated by the blue arrows and is used to increase or decrease the cross section of the duct geometry. One design variable per node is needed for this morph shape which includes four or eight control points with radial vectors that allow for the change in the cross section. The total number of morph shapes and thus the number of design variables for shape optimization is 66.

Based on the automatically computed morph shapes, shape optimization is performed further minimizing the power dissipation, with the same constraints on the mass flow rates at the exits and constraint on the volume of the design to be constant with a small tolerance. The initial and optimal geometries during this stage of optimization are shown in Figure 12. The distribution of the magnitude of the grid displacement that corresponds to the optimal geometry is also shown. The mass balance constraint at the exits is satisfied at this stage with a tolerance of 1%.
Figure 12. Shape optimization: initial and optimal shape


The optimal design of the HVAC with streamlines and velocity vectors is shown in Figure 13, together with the continuing velocity vectors inside the cabin.
Figure 13. Optimal design


References

[1] Sandboge, R., Papadimitriou, D., and Reddy, M., “Shape and Topology Optimization in Computational Fluid Dynamics Including Heat Transfer Using Gaussian Processes and Adjoint Methods,” 2021. https://doi.org/10.2514/6.2021-1892.

[2] Svanberg, K., “MMA and GCMMA-two methods for nonlinear optimization,” vol, Vol. 1, 2007, pp. 1–15.