Industrial application of topology optimization for combined conductive and convective heat transfer problems

Abstract

This paper presents an industrial application of topology optimization for combined conductive and convective heat transfer problems. The solution is based on a synergy of computer aided design and engineering software tools from Dassault Systèmes. The considered physical problem of steady-state heat transfer under convection is simulated using SIMULIA-Abaqus. A corresponding topology optimization feature is provided by SIMULIA-Tosca. By following a standard workflow of design optimization, the proposed solution is able to accommodate practical design scenarios and results in efficient conceptual design proposals. Several design examples with verification results are presented to demonstrate the applicability.

1 Introduction

Topology optimization is a generative design tool for conceptual structural design, pushing the limit of product performance based on computer simulation and optimization technologies. It originates from the area of mechanical structural design (Bendsøe and Kikuchi 1988; Bendsøe and Sigmund 2003) (e.g. aircraft, automotive, etc.) where a lightweight structural layout is desired that satisfies static or dynamic design requirements. It now has been broadened to industrial application for multi-disciplinary problems.

This paper presents an industrial solution to design topologically optimized structures for combined conductive and convective heat transfer problems. The solution is a synergy of various technologies, including Computer Aided Design (CAD), Finite Element Analysis (FEA) and topology optimization, which are enabled by Dassault Systèmes softwares Catia, SIMULIA-Abaqus (2015) and SIMULIA-Tosca (2015), respectively. The proposed solution can be applied to the industrial design of electronic devices, heating appliances, combustion engines and electric motors etc., when thermal management becomes critical besides other design requirements such as light weight, high stiffness and strength.

Fig. 1

Topology optimization with a constant and a design-dependent convection: a design domain and boundary conditions; b optimized density distribution using a constant and design-independent convection; c optimized density distribution using a design-dependent convection

The considered physical problem herein is steady-state heat transfer under convection where thermal energy is transferred through (i) conduction inside the structural body and (ii) convection at the structure-fluid interface. For topology optimization of such a problem, one key issue is to characterize the convection load in the optimization process. With the density-based topology optimization approach, the convection surface is loosely defined at the beginning and changes until the overall topology has been settled. Since convection is a surface controlled phenomenon, it is important to find a way to effectively interpolate the convection into the design domain. In previous works, various ways have been proposed to model the convection: Sigmund (2001) uses a simplified convection model where side convection is indirectly modeled by conduction into void regions; Yin and Ananthasuresh (2002) use a density-based peak interpolation function; Yoon and Kim (2005) introduce a special type of parameterized connectivity between elements; Bruns (2007) suggests to interpolate the convection coefficient (also known as film coefficient) as a function of the density-gradient; Iga et al. (2009) use a density-based smeared-out Hat-function, which also try to take variation in the strength of the convective heat transfer into account; Dede et al. (2015) adopt a solution similar to Iga et al. (2009) to design and manufacture heat sinks subject to jet impingement cooling including spatial variability of the convection coefficient. Another approach is the level set method based optimization. In the works by Ahn and Cho (2010) and Coffin and Maute (2015), the convection is defined precisely at the structure-fluid interface.

The topology optimization feature for combined conductive and convective heat transfer problems in SIMULIA-Tosca follows the idea which was initially suggested by Bruns (2007) and investigated by Alexandersen (2011a, (Alexandersen 2011b)). A simplified engineering approach is adopted by assuming a design-dependent film coefficient in the optimization process, which is readily compatible with the thermal analysis in SIMULIA-Abaqus. It should be mentioned that it is possible to apply topology optimization to the full conjugate heat transfer problem (see (Yoon 2010) for forced convection and (Alexandersen et al. 2014) for natural convection). However, this increases the computational cost significantly and is currently not desirable in industrial settings.

The remaining content of this paper is organized as follows. In Section 2, the governing equation and finite element model of steady-state heat transfer under convection are first introduced. Then, the importance of applying a design-dependent convection in topology optimization is highlighted using a demonstration example. In addition, an industrial workflow of design optimization using CAD and CAE softwares is stated. In Section 3, several industrial design examples are presented and discussed with verification results. Conclusions are stated in Section 4.

2 Problem formulation

2.1 Steady-state heat transfer under convection

The governing equation for the computational domain O is the steady-state heat transfer equation:

$$ -\frac{\partial}{\partial x_{i}} \left( k \frac{\partial T}{\partial x_{i}} \right) = 0 \text{ in } \Omega, $$

(1)

where T is the temperature field, x _i are the spatial coordinates and k is the thermal conductivity. The boundary is split into disjoint subsets ${\Gamma } = {\Gamma }_{\text {flux}} \cup {\Gamma }_{\text {ins}} \cup {\Gamma }_{\text {conv}} \cup {\Gamma }_{\text {temp}}$, on which the following boundary conditions are imposed:

$$\begin{array}{@{}rcl@{}} -k \frac{\partial T}{\partial x_{i}}n_{i} = q_{0}\hspace*{3.2pc} & & \qquad\qquad\text{on } {\Gamma}_{\text{flux}} , \end{array} $$

(2)

$$\begin{array}{@{}rcl@{}} -k \frac{\partial T}{\partial x_{i}}n_{i} = 0\hspace*{3.5pc} & & \qquad\qquad\text{on } {\Gamma}_{\text{ins}} , \end{array} $$

(3)

$$\begin{array}{@{}rcl@{}} -k \frac{\partial T}{\partial x_{i}}n_{i} = h(T-T_{\text{ref}}) & & \qquad\qquad\text{on } {\Gamma}_{\text{conv}}, \end{array} $$

(4)

$$\begin{array}{@{}rcl@{}} T = T_{0}\hspace*{3.1pc} & & \qquad\qquad\text{on } {\Gamma}_{\text{temp}}, \end{array} $$

(5)

where (2) is the prescribed surface flux condition, (3) is the insulation condition, (4) is the boundary convection condition, and (5) is the prescribed temperature condition. Furthermore, q ₀ is the prescribed surface flux, h is the convection coefficient, T _ref is the reference temperature of the ambient fluid and T ₀ is the prescribed temperature.

By discretizing the governing equation using finite elements, a system of linear equations is obtained as follows:

$$ \left( \mathbf{K} + \mathbf{H} \right) \mathbf{t} = \mathbf{f} + \mathbf{f}_{h}, $$

(6)

where K is the conductivity matrix, H is the convection matrix, t is the vector of nodal temperatures, f is the flux vector arising from (2) and f _h is the convection vector arising from (4).

2.2 Design-dependent convection in topology optimization

Heat transfer by convection happens when solid objects are in contact with a fluid and the energy is transferred to or from a surrounding fluid across the solid-fluid interface due to a temperature difference. Generally, the efficiency of heat transfer is determined by various factors, including the property and speed of the moving fluid and the size, shape and properties of the solid object. It is possible to simulate a conjugate heat transfer process using thermofluidic modeling, which can also be leveraged in topology optimization (Yoon 2010; Alexandersen et al. 2014). However, it requires intensive computational resources and it is currently not desirable for industrial applications.

Another simplified way to characterize the heat transfer by convection is to define an average convection heat transfer coefficient and a convective flux at the interface by (4). The convection coefficient can be estimated or found in engineering tables for different solid-fluid interactions. In the early design phase when designers neither know the exact working environment nor the geometry of the component, using an approximation for heat-transfer is common engineering practice. Besides, fast turnaround time in the early conceptual design phase is crucial, making a simplified simulation approach acceptable at the cost of lack in accuracy. Therefore, a topology optimization approach considering the simplified approach is of industrial importance.

A topology optimization feature for thermal conductive and convective problems is implemented in SIMULIA-Tosca. It follows the ideas proposed by Bruns (2007) and later studied by Alexandersen (2011a, (Alexandersen 2011b)): that both the conductivity and the convection coefficient are assumed design-dependent and interpolated during the optimization process. In the current implementation, the conductivity k _i for the element i is defined as k _i = g _k (? _i )k ₀, where g _k is the interpolation function for thermal conductivity, ? _i is the element density and k ₀ is the specified material conductivity. The convection coefficient h _{i j} at the interface between element i and j is defined based on the density differences between adjacent elements. Mathematically, it is expressed as h _{i j} = g _h (? _i -? _j )h ₀, where g _h is the interpolation function for the convection coefficient and h ₀ is the specified ambient convection coefficient. However, due to confidentiality, the actual functions g _k and g _h are not given in this paper. Interested readers are referred to a similar work by Alexandersen (2011a) for more technical details on using interpolation strategies for thermal problems.

The importance of using design-dependent convection in topology optimization is verified using a demonstration example. As shown in Fig. 1a, a rectangular block of dimension 250×250×40 mm is setup as the design domain. It is subject to a prescribed temperature T=100 ^°C at the central region of the left surface (red color) and convection based on a reference ambient temperature of $T_{\text {ref}}=20{~}^{\circ }\mathrm {C}$. The conductivity of the solid material is assumed to be k=0.385 W/(mm·K) and the convection coefficient is h=10^-5 W/(mm²·K). The model is discretized with 100×100×4 linear hexahedral elements. The optimization objective is to maximize the reaction flux¹ at the region where prescribed temperature applies. In addition, the problem is subject to a volume constraint, which must be less or equal to 40 % of the design domain. Figure 1b shows the optimized density distribution² by using a design-independent convection, where a constant convection coefficient is assumed over the outer boundary of the design domain during the optimization. The corresponding validation result in Table 1 shows that it achieves a total reaction flux of 138.5 W at the prescribed-temperature region. As a comparison, the optimized design with a design-dependent convection is given in Fig. 1c. It not only shows different topology and shape from the previous design but also possesses a larger surface area and a better thermal conductive behavior with a total reaction flux of 144.2 W at the prescribed-temperature region. Note that both designs are obtained under the same optimization parameters, except for the convection scheme. Thus, it can be concluded that even for the present simple heat transfer problem, by using the design-dependent convection scheme has a significant impact on both the optimized structural layout as well as the objective value.

Table 1 Verification of the optimized designs obtained by different convection schemes

2.3 An integrated CAD and CAE workflow

The proposed solution is based on an integrated Computer Aided Design (CAD) and Engineering (CAE) workflow enabled by Dassault Systèmes. As shown in the flow chart in Fig. 2, a geometric model is initially prepared in a CAD software, such as Catia or Solidworks, and then a finite element model is created in SIMULIA-Abaqus for topology optimization. Upon setting up the optimization goals and constraints with thermal related design responses, a simulation based topological design process is executed using SIMULIA-Tosca. For each design iteration, the considered problem of steady-state heat transfer under convection is simulated using SIMULIA-Abaqus. After that, the optimization module performs an adjoint sensitivity analysis and uses a gradient-based optimizer to search for an optimized design. Once the optimization process converges by satisfying certain stopping criteria, the optimized density model is post-processed by iso-surface extraction and smoothing and then transferred into a surface model. The post-processing operations may slightly change the geometry as well as the structural performance of the design. However, this is an indispensable step in practice for any downstream applications, which require a CAD model rather than a density model as input. The post-processing step and its impact are discussed in the numerical example sections. In addition, a shape conforming finite element model is generated and tested for performance verifications.

Fig. 2

Workflow of topological designing using Dassault Systèmes software solutions

3 Design examples

The present section shows several examples demonstrating the applicability of the proposed solution. Manufacturing restrictions such as minimum member size, symmetry, casting and extrusion constraints are included in the topology optimization process. These constraints are enabled by SIMULIA-Tosca (2015) and they are similar to well-established methods as proposed by Guest et al. (2004), Zhou et al. (2015) and Gersborg and Andreasen (2011). The minimum-member-size constraint is imposed to control the size of the thinnest structural member. The symmetry constraint not only leads to a symmetric structural design which is aesthetically desirable but can also compensate numerical errors in optimization when the problem is subject to symmetric design domain and boundary conditions but asymmetric meshing (c.f. Section 3.1). The casting and extrusion constraints are applied to avoid undercut along the pulling direction for casting and to prevent possible nucleation of internal holes. The latter will cause unrealistic energy loss based on the simplified engineering convection scheme considered. Note that the optimized structure may not be castable if more than one casting constraint is used. Some key results and statistics of the workflow, including the optimized density distributions, the convergence curves, the corresponding smoothed designs and verification results are provided.

3.1 Topology optimization with different types of finite elements

One prerequisite for general industrial application of topology optimization is to support different types of finite elements. The present section demonstrates the robustness of SIMULIA-Tosca on topologically designing a heat conductor using different types of elements. The initial geometry of the model and boundary conditions are given in Fig. 3a, where a cube is subject to a heat flux of magnitude 2.0 W/mm ² at the central region (purple color) of the bottom surface, a prescribed temperature $T_{0}=0{~}^{\circ }\mathrm {C}$ at different locations (red color) on the top surface and surface convection with a reference temperature of $T_{\text {{ref}}}=20{~}^{\circ }\mathrm {C}$ and h=10^-5 W/(mm²·K) over the design domain (brown color). The dimensions of the cube are 40×40×40 mm and the conductivity is k=0.237 W/(mm·K). The design problem is to minimize the thermal compliance subject to a volume constraint of 20 % of the design domain. The thermal compliance is defined as a summation of the product of the equivalent nodal heat flux and temperature values at the external flux region. The convergence criterion is a measure based upon (i) the change in the objective function to be less than or equal to 0.1 % (ii) the change in the density to be less than 0.1 % between a current and a previous optimization iteration.

Fig. 3

Topology optimization using different types of elements: a boundary conditions; b–q (for each row from left to right) element type, optimized density distribution (refer to Fig. 1 for the color scale) and smoothed designs (snapshots from different views). Element types (from top to bottom) 4-node tetrahedron (DC3D4), 10-node tetrahedron (DC3D10), 8-node linear hexahedron (DC3D8) and 20-node linear hexahedron (DC3D20), respectively

The cube is discretized with four different types of finite elements, namely linear tetrahedrons (DC3D4), quadratic tetrahedrons (DC3D10), linear hexahedrons (DC3D8) and quadratic hexahedrons (DC3D20), respectively. The optimized designs are compared in Fig. 3b–q. Each row of figures shown from left to right lists the element type, the density distribution and the smoothed designs from isometric and top views, respectively. All of the four optimizations converge to the same topology and almost the same shape. The statistics of the optimization, including the number of elements and nodes per mesh, and verification results of the four smoothed designs are listed in Table 2. The verifications are based on remeshing the smooth designs with quadratic tetrahedron (DC3D10) meshes. All of the four designs have an average temperature around 80 ^°C at the external flux region. The slight difference in the volume is due to postprocessing, e.g. the iso-surface cut from the density model and smoothing. The performances of the designs obtained by the tetrahedron meshes are better than those by the hexahedron meshes because the tetrahedron meshes allow for more detailed designs with more elements. The convergence curves of all designs are shown in Fig. 4. Note that two additional symmetry constraints are imposed for the tetrahedron meshes in order to ensure symmetric results as shown in Fig. 3c–e, g–i. Otherwise the optimized results are slightly asymmetric due to asymmetric meshing of the design domain.

Fig. 4

Optimization iteration history using four different element types as shown in Fig. 3

Table 2 Statistics of the optimization in Fig. 3 and verifications (using DC3D10 meshes) for different smoothed designs

The smoothed design proposals shown in Fig. 3 exhibit a structural topology connecting the heat source to the four sinks (T=0 ^°C) at the middle regions of the edges on the top surface. The other four heat sinks at the corners as shown in Fig. 3a have a longer distance from the heat source and thus, they are less efficient in energy transfer. Verification shows that only 2 % of the energy is transferred through convection in this example. Hence, the problem here is conduction-dominated. Nonetheless, it demonstrates the flexibility of SIMULIA-Tosca in handling different types of finite elements and yielding consistent results irrespective of continuum element types for heat-transfer based topology optimization.

3.2 Topological design of an electric motor cover

3.2.1 Optimized design with thermal design responses

In this section, the proposed solution is applied to design an electric motor cover. The model is sketched in Fig. 5, where a heat producing motor is surrounded by a solid and non-designable casing illustrated by red in Fig. 5a and a uniform heat flux of magnitude 0.02 W/mm ² is generated by the motor illustrated by orange in Fig. 5b. The non-design domain comprises two horizontal mounting brackets as well as bearings for the electric rotor. The remaining domain constitutes the design domain shown in green. The design goal is to conceive a cover structure which transfers the heat efficiently to the surroundings via conduction and convection. This model has a length of 170 mm and a radius of 70 mm. It is discretized with a mesh consisting of 46574 linear hexahedral elements, 2522 linear wedge elements and 1474 linear tetrahedral elements. Other parameters are set as: k=0.385 W/(mm·K), h=2×10^-5 W/(mm²·K) and $T_{\text {ref}}=20{~}^{\circ }\mathrm {C}$. The optimization problem is to minimize the thermal compliance subject to a volume constraint of 50 % of the entire structure, including both design and non-design domains.

Fig. 5

Geometric model and boundary conditions of an E-motor cover: a design domain (green color) and non-design domain (red color); b external heat flux applied at the inner wall and convection defined over the structural boundary

By applying casting constraints along the radial direction of the motor’s central axis, two conceptual design proposals with different minimum member sizes are presented in Figs. 6 and 7, respectively. For comparison, the optimized density distribution, the corresponding smoothed design, the verification results showing the temperature distribution as well as the optimization history curves are given for both designs. The first design as shown in Fig. 6 is obtained by applying a minimum-member-size restriction of 12 mm. For such a design, spike-shape structures appear over the electric motor cover in order to effectively transfer the heat to the environment. By enlarging the minimum length scale to 20 mm, thicker wall structures are observed in the optimized design as shown in Fig. 7. The detailed verification results for average temperature at the flux region as well as the maximum and the minimum temperature over the structure are recorded in Table 3. The verification reveals as expected that the design with a smaller member size achieves a better heat transferring capability compared to the other constrained design, though the volume ratios of two verification models are slightly different.

Fig. 6

Topological design I of an E-motor cover for a small minimum member size: a density distribution; b smoothed design; c temperature distribution (^°C); d optimization history

Fig. 7

Topological design II of an E-motor cover for a large minimum member size: a density distribution; b smoothed design; c temperature distribution (^°C); d optimization history

Table 3 Validation value of E-motor cover designs

3.2.2 Optimized designs with different material and convection coefficients

For the combined conductive and convective heat transfer problems, the dimensionless Biot number $B=\frac {h\cdot L_{c}}{k}$, where L _c is a reference length, quantifies the ratio of heat transfer resistances at the surface and inside of a structural body. For topology optimization of such a problem, different Biot numbers will result in different optimized structural layouts (Alexandersen2011a, (Alexandersen 2011b); Coffin and Maute2015). Generally, a low B value leads to a design with pronounced structural features, such as elongated fins, for which the surface area is maximized to obtain an efficient heat dissipation through convection. In contrast, a high B value results in an optimized design with shorter and thicker fins since conductive resistance now becomes the bottleneck and such resistance through the fins must be minimized. Such phenomena are reproduced here by designing the E-motor cover with copper and steel using different convection strengths. Figure 8 compares two design proposals, where Fig. 8a is the same as that in Fig. 6 with k=0.385 W/(mm·K), h=2×10^-5 W/(mm²·K) and Fig. 8b is obtained with a lower conductivity k=0.065 W/(mm·K) and a higher convection coefficient h=5×10^-3 W/(mm²·K). For the latter case with a high B value, the conductive material accumulates around the heat source of the E-motor cover in layers of fins along its central axis. Note that both optimization processes are executed for the same volume fraction and with the same minimum member size. The length of the E-motor cover is used as the reference length L _c =170 mm for calculating the Biot numbers. The optimization results here are in an agreement with the previous studies in (Alexandersen2011a, (Alexandersen 2011b); Coffin and Maute2015). Readers are referred to the references for further discussions on the impact of applying different materials and convection coefficients in topology optimization.

Fig. 8

Topological design of E-motor cover with different ratios of convection coefficient and thermal conductivity: a B=0.9×10^-3 and b B=13.1

3.3 Topological design of heat sinks

3.3.1 Optimization with thermal and dynamic design responses

Heat sinks are passive heat exchangers that are used widely in electronic devices to cool key components such as central processing units (CPUs) or graphics cards in a computer. In this section, several heat-sink designs are presented. Their thermal behaviors and fundamental structural frequencies are compared. The initial geometry and boundary conditions are given in Fig. 9. It consists of a design domain (the upper cubic part) of dimension 100×100×50 mm and a non-design domain (the bottom base) of dimension 40×40×10 mm for assembly purpose. A uniform surface heat flux of magnitude 0.1 W/mm² is applied at the bottom surface and a forced convection exists over the structural boundaries. The optimization goal is to minimize the thermal compliance s.t. a volume constraint of 30 % of the design domain. The model is discretized into 264192 linear hexahedral elements for finite element analysis. For the solid material, the Young’s modulus of 70 GPa, the density of 2.8×10^-6 kg/mm³ and the Possion ratio of 0.35 are assumed. Other parameters are set as follows: k=0.237 W/(mm·K), h=1×10^-4 W/(mm²·K) and $T_{\text {ref}}=20{~}^{\circ }\mathrm {C}$.

Fig. 9

Topology optimization of a passive heat sink: design domain (brown color), non-design domain (blue color) and boundary conditions

Figure 10 illustrates a reference heat-sink design and the corresponding simulation results. The sink possesses a desired mass fraction of 30 % as specified. It has a fundamental structural frequency of 50.3 Hz, an average temperature of 61.4 ^°C at the bottom surface, a maximum temperature of 62.2 ^°C and a minimum temperature of 44.1 ^°C. For comparison, three other design proposals with the same volume fraction and a similar length scale are obtained using topology optimization. The thermal behavior and the fundamental structural frequency of each heat-sink design are recorded in Table 5, including the average temperature of the bottom surface (where the heat flux applies), the maximum and minimum temperature of the overall structure as well as the fundamental structural frequency.

Fig. 10

A reference heat-sink design: a geometric model; b,c temperature distribution (^°C); d mode shape of 1^st eigenfrequency

Figure 11 shows the first optimized heat sink (Design I) and its verification results. This design is obtained by applying three manufacturing constraints (extrusion, casting and symmetry) in the topology optimization process. The direction of the extrusion is out of plane and the casting direction is in the 45-degree direction on each half of the structure. The overall design process converges in 26 design iterations with an active volume constraint. The second design proposal (Design II) is given in Fig. 12, which is obtained by applying one symmetry constraint and one casting constraint in the upward direction normal to the base. This design possesses spikes and thin-walled structures which follow the casting direction and spread over the base. The optimization process takes 28 design iterations to converge with an active volume constraint.

Fig. 11

Heat-sink design I: a smoothed design; b,c validation of temperature distribution (^°C); d mode shape of 1st eigenfrequency

Fig. 12

Heat-sink design II: a smoothed design; b,c validation of temperature distribution (^°C); d mode shape of 1st eigenfrequency

To evaluate the performance difference between the density model and the final CAD design, the optimization valueof Designs I and II from SIMULIA-Tosca are recorded in Table 4. Note that the minimum temperature in the design domain appears in the void region, which equals to the ambient reference temperature. Comparing to Table 5, the verification results slightly differ from the optimization value as shown in Table 4 for several reasons. First, design-dependent convection coefficients are assumed inside the design domain in optimization, while in verification a constant convection coefficient is considered over the structural surface of a solid-void CAD model. Second, the verification model (CAD model) deviates from the density model after the iso-surface extraction, smoothing and local surface modifications. The differences in the volume fraction between the verification and optimization models are acceptable from an engineering point of view. In this example, Design II exhibits a better heat transfer capability as its average temperature at the flux regions is lower than that of Design I. Besides, as a by-product of this optimization, the maximum temperature of Design II is also lower than that of Design I.

Table 4 Optimization value of heat-sink designs (Density model)

Table 5 Verification of the smoothed heat-sink designs using a simplified engineering model

In practice, it is common to have both static and thermal design requirements for the same design. Here, a third design proposal is obtained by maximizing the fundamental structural frequency³ of the heat sink subject to structural volume and thermal constraints. The constraints require that the maximum temperature of the heat sink must not exceed 60.0 ^°C and the final structural volume is less or equal to 30 % of the initial geometry. Besides, the same manufacturing restrictions as those in Design I are applied. Maximizing the fundamental frequency ensures a generally stiff design which does not require definition of specific static load cases but may also be applied to ensure that the fundamental frequency stays above the rotational frequency of e.g. nearby cooling fans. Figure 13 shows the smoothed design, the corresponding temperature distribution and the 1^st eigenmode shape. Comparing to Design I, this heat sink possesses shorter fins with major material accumulating around the base. It achieves a fundamental frequency of 147.9 Hz and a maximum temperature of 73.6 ^°C in the post-processed design. Note that, in the optimization process the temperature constraint $T_{max}\leq 60{~}^{\circ }\mathrm {C}$ is strictly satisfied. However, the final temperature of the smooth result deviates from the constrained value due to postprocessing the optimized density model (with intermediate densities) and transferring it into a solid-void CAD model.

Fig. 13

Heat-sink design III: a optimized design; b,c temperature distribution (^°C); d mode shape of 1st eigenfrequency

3.3.2 Verification with a thermofluidic model

The reference heat sink and optimized Design I in Section 3.3.1 are further verified⁴ with thermofluid based conjugate heat transfer in COMSOL-Multiphysics (2015). As shown in Fig. 14a, the heat sink is placed at the bottom of a computational domain, where a laminar air inflow with an average velocity 3.9 m/s and initial temperature 20 ^°C is applied from the left-front inlet to the right-back open boundary. The other four sides of the domain are thermally insulated and no fluid penetrates. The material property of the heat sink and the external heat flux are the same as those in the previous verification using SIMULIA-Abaqus.

Fig. 14

Verification of the reference heat sink using thermofluidic modeling

In this thermofluidic verification, the average inlet fluid velocity is tuned such that the average heat exchange efficiency due to the structure-fluid interaction is equivalent to that in the optimization. Specifically, an effective convection coefficient is defined based on the following equation:

$$ h^{\prime}=\frac{1}{A}{\int}_{\!\!\!\!{\Gamma}}{\frac{q}{(T-T_{\text{ref}})}d{\Gamma}}, $$

(7)

where A is the total surface area of the structure-fluid boundary G subject to heat convection(excluding the bottom surface), T is the temperature, T _ref is the ambient reference temperature and q is the flux from the structure to the fluid. As recorded in Table 6, with the fluid velocity 3.9 m/s, the effective convection coefficient $h^{\prime }$ for the heat-sink Design I is 1.00×10^-4 W/(mm²·K), which is the same as the convection coefficient value used in the optimization.

Table 6 Verification of the smoothed heat-sink designs using thermofluidic modeling, inlet fluid velocity 3.9m/s

The steady-state fluid streamlines and temperature distribution over the reference heat sink are given in Fig. 14a–b, respectively. As a comparison, Fig. 15a–b shows the corresponding simulation results of the optimized design. Comparing to the verification results as shown in Figs. 10 and 11, which use a simplified convection model and a uniform ambient temperature, the temperature distribution over the heat sink here increases gradually in the direction of the fluid flow.

Fig. 15

Verification of heat-sink Design I using thermofluidic modeling

The thermal behaviors of two designs, including the average temperature at the flux region, the maximum and minimum temperature over the heat sink, the effective convection coefficient under the fluid of average velocity 3.9m/s as well as the average pressure drop between the inlet and outlet boundaries, are evaluated and listed in Table 6. Despite the temperature values being different from those in Table 5, the optimized Design I still exhibits a better thermal performance here than the reference design by using a thermofluidic model. The former has a lower average temperature of 63.5 ^°C at the flux regions comparing to 83.5 ^°C. Besides, the maximum and minimum temperatures of the Design I are 20.0 and 18.2 ^°C lower than the other, respectively. Another important metric is the pressure drop of the fluid between the inlet and outlet, which relates to pumping power. The verification results in Table 6 show that the fluid passing the reference design has a lower average pressure drop of 0.6 Pa between the inlet and the outlet boundaries, compared to 0.8 Pa for the optimized Design I. This increased pressure drop can probably not be avoided since better cooling requires larger surface area and hence causes larger pressure drop. Anyway, the considered simplified convection model cannot take this metric into account and hence further studies on such effects require a full-blown thermofluidic analysis model.

To this end, both verifications of using either a simplified engineering model or an advanced thermofluidic one indicate that the proposed design optimization solution is able to yield efficient structural designs with optimized thermal behaviour for combined conductive and convective heat transfer problems.

3.3.3 Design with different discretizations

This section studies the impact of the mesh discretization on the optimization results. The heat-sink Design I in the previous section is used as the example for comparison. As shown in Fig. 16, three sets of optimized designs are obtained by using finite element meshes of 64000, 264192 (same as Design I in Section 3.3.1) and 560000 hexahedral elements, respectively. The thermal behaviours of the smoothed designs are evaluated and recorded in Table 7. The design based on the coarse mesh exhibits the highest average temperature 60.3 ^°C at the flux region among the three. Some geometric details, such as the fins at the bottom corners of the design domain, are missing comparing to the others due to insufficient design freedoms. The normal and fine meshes yield similar optimized geometries, except that the latter has a smoother boundary and slightly different finger sizes. However, the design by the fine mesh has a higher average temperature (57.9 ^°C) at the flux region than that of the normal one (56.6 ^°C). It is because that the design by the normal mesh, as shown in Fig. 16b, has intermediate densities at the finger tips and the minimum length scale is not strictly fulfilled there. As a result, it leads to sharp and slender finger tips in the smoothed design as shown in Fig. 16e, which allows for a more efficient heat dissipation.

Fig. 16

Heat-sink design with different discretizations: (from the left column to right) meshes of 64000, 264192 (same to Design I in Section 3.3.1) and 560000 hexahedral elements; (top row) optimized density models (refer to Fig. 1 for the color scale); (bottom row) smooth designs

Table 7 Verification of the smoothed heat-sink designs in Fig. 16

4 Conclusions

An industrial solution for topology optimization of combined conductive and convective heat transfer problems is implemented in Dassault Systemes software solutions and numerical industrial experiments show promising results. The topology optimization feature is provided by SIMULIA.Tosca, where a simplified engineering model is utilized to model the design-dependent convection during the optimization process. Such a strategy allows for a seamless integration with the static-state heat simulation module in SIMULIA.Abaqus. Several design examples with different problem setups are given to demonstrate its applicability. Verification results of each optimized design are derived by either steady-state heat transfer using the engineering convection scheme or conjugate heat transfer with a thermofluidic model. Both studies indicate that the proposed solution is able to yield efficient design proposals for various practical design scenarios.