# On the Topology Optimization of Elastic Supporting Structures under Thermomechanical Loads

Authors: Jie Hou; Ji-Hong Zhu; Qing Li
Journal: International Journal of Aerospace Engineering
 (2016)

---

## Abstract

This paper is to present a thermomechanical topology optimization formulation. By designing structures that support specific nondesignable domain, optimization is to suppress the stress level in the nondesignable domain and maintain global stiffness simultaneously. A global stress measure based onp-norm function is then utilized to reduce the number of stress constraints in topology optimization. Sensitivity analysis employs adjoint method to derive the global stress measure with respect to the topological pseudodensity variables. Some particular behaviors in thermomechanical topology optimization of elastic supports, such as the influence of different thermomechanical loads and the existence of intermediate material, are also analyzed numerically. Finally, examples of elastic supports on a cantilever beam and a nozzle flap under different thermomechanical loads are tested with reasonable optimized design obtained.


---

## Body

## 1. Introduction

Optimization design of thermoelastic structures is an important problem in aeronautics and aerospace products such as turbine engine components and thermal protection systems (TPS). As one of the most challenging topics in topology optimization, many difficulties are to be settled before effective engineering application.One the one hand, literatures dealing with structural configurations design for different heat conduction purposes have been present using topology optimization, where the total potential energy was usually minimized as the design objective. For example, Li et al. [1] optimized the heat conduction path using the typical evolutionary structural topology optimization. The ineffective materials were gradually removed from the design domain iteratively. The final designs turned out to be some tree-like structures for volume-to-point heart conduction. Later, density based topology optimization methods were used to obtain similar optimized designs [2, 3]. Gao et al. [4] recognized the heat sources as design-dependent thermal loads; some more effective and clear optimized designs were finally obtained. Later, the design-dependent effect of heat conduction and convection in homogenization based topology optimization was further investigated [5]. Recently, Dirker and Meyer [6] extended the heat conduction topology optimization to the internal cooling system design.One the other hand, topology optimization of coupled thermoelastic structures was generally considered as mean compliance minimization problems with different solutions in earlier works (e.g., [7–13]). Gao and Zhang [14] then extended such formulation to thermoelastic structures of multiphase materials. The concept of thermal stress coefficient (TSC) defined by Young’s modulus and thermal expansion coefficients was proposed to characterize the dependence of thermal stress and thermal expansion coefficient upon element pseudodensities. Deaton and Grandhi [15] investigated the compliance design of the supporting structures in thermoelastic environment. Nonintuitiveness in design space and the design dependency that occurs with thermal loading were presented and different topology optimization formulations were used to solve the problem. A recent literature survey by Deaton and Grandhi [16] has summarized the advances and applications of thermoelastic topology optimization.Actually, eligible design of thermoelastic structures should meet both requirements of structural stiffness and stress to prevent structural failure. As a result, recent efforts have been continuously devoted to highlighting thermoelastic problems with the consideration of compliance, strain energy, and maximum stress. P. Pedersen and N. L. Pedersen [17, 18] proposed strength optimization minimizing the maximum Von-Mises stress for thermoelastic structures, where thermal stress was generated by thermal expansion of clamped structures. Zhang et al. [19] presented a topology optimization formulation of elastic supports for thermoelastic structures and investigated the strength design by directly minimizing the maximum Von-Mises stress in the nondesignable subregion. They further analyzed the differences of elastic strain energy and mean compliance of thermoelastic structures. Recently, Deaton and Grandhi [20] presented a topology optimization method with combined mechanical and thermal loads. A relaxation technique and a modified p-norm function were utilized to remove the singularity phenomenon and aggregate the large number of stress constraints.In this paper, we continue to study the thermoelastic structures topology optimization suppressing stress level and maintaining global stiffness simultaneously. As shown in Figure1, a typical thermoelastic continuum structure undergoes coupled thermal and mechanical loads. Subregion Ω is nondesignable and Ω
s is the elastic supporting structure assigned as the design domain. Large stress level appears in the nondesignable domain due to the coupled effect of thermal expansion and mechanical loads. Design objective here is to find proper structural configuration of the elastic supports that suppress the stress in the supported nondesignable domain and simultaneously maintain global stiffness evaluated by structural compliance. As the stress constraints are only applied on the nondesignable domain, the singularity problem in stress-based topology optimization is actually not involved here.Figure 1
Topology optimization of thermoelastic structure with elastic support.To have an in-depth understanding of the thermoelastic problems with elastic supports, we perform a detailed analysis of a three-bar truss analytical model to present some illustrative phenomena during optimization. The behaviors of structural compliance and stress are discussed, respectively. Numerical examples of elastic supports on a 2D cantilever beam and 3D nozzle flap are finally presented to verify the validity of the presented formulation.


## 2. Three-Bar Truss Analytical Model

The investigation on a truss-frame structure can reveal the underlying scheme of continuum topology optimization (e.g., [21, 22]). In this section, the analytical solutions of a three-bar truss system are derived under different thermomechanical loads with regard to the stress constraints. The structural compliance of thermoelastic problems is evaluated with further details discussed here.As sketched in Figure2, the truss system is comprised of three equal-length bars, that is, E
(
1
), E
2, and E
3. The bars’ cross-sectional areas are assigned as x
1, x
2, and x
3, respectively. One end of those bars is fixed and the other ends are hinged together to node 4. The mechanical loads applied on node 4 can be decomposed into two components f
m
x and f
m
y. Meanwhile, all three bars undergo a temperature increase Δ
T and the corresponding thermal expansion will cause thermal stresses. E
1 and E
2 constitute the elastic support and the cross-sectional areas, that is, x
1 and x
2, which are assigned as the design variables. E
(
3
) is chosen as the nondesignable domain which means x
3 is fixed. In this example, the continuous design variables are used to evaluate the behaviors of the structural response.Figure 2
Three-bar truss system with elastic support.According to the thermoelastic theory, the thermal load in each bar can be calculated as(1)
F
i
t
h
=
E
α
Δ
T
x
i
,where E is Young’s modulus, α is the thermal expansion coefficient, Δ
T is the temperature increase, and x
i is the cross-sectional area of the corresponding bar.The strain and stress of the barE
i can be expressed as(2)
ε
i
=
U
i
x
¯
L
,
σ
i
=
E
ε
i
m
=
E
ε
i
-
ε
i
t
h
=
E
U
i
x
¯
L
-
α
Δ
T
,where U
i
x
¯ is the displacement with respect to the elemental x-axis and ε
i
m and ε
i
t
h are corresponding mechanical and thermal strain.The elastic strain energyϕ is defined as the potential mechanical energy in the elastic body, which is written as(3)
ϕ
=
1
2
∫
ε
m
T
D
ε
m
d
V
=
1
2
∫
ε
-
ε
t
h
T
D
ε
-
ε
t
h
d
V
=
1
2
∫
ε
T
D
ε
d
V
-
∫
ε
T
D
ε
t
h
d
V
+
1
2
∫
ε
t
h
T
D
ε
t
h
d
V
,where D is the elastic matrix and ε is the total strain vector consisting of mechanical and thermal items ε
m and ε
t
h.Then, we consider the mean compliance of structure(4)
C
=
1
2
F
T
U
=
1
2
F
T
K
-
1
F
,where F and U are nodal load and displacement vectors. K is global stiffness matrix.According to the existing works, the strain energy minimization is more beneficial for the stress reduction while the mean compliance reflects the structural overall stiffness [19]. The purpose of this study is to find the optimized configuration of elastic supports preventing large deformations resulting from thermal and mechanical load. Naturally, the mean compliance minimization is selected as the design objective throughout this paper. Besides, two design constraints, that is, volume constraint and stress constraint, are involved in this study. The formulation of the optimization can be expressed as(5)
Find:
x
1
,
x
2
Min:
C
s.t.:
x
1
+
x
2
L
≤
V
¯
;
σ
3
≤
σ
¯
.To demonstrate the relation between the design variables and the global compliance, parameters listed in Table 1 are used.Table 1
Constant parameters list.
Young’s modulusE (Pa)
2
×
1
0
11
Thermal expansion coefficientα (°C−1)
1
×
1
0
-
5
Length of barL (m)
1
Cross-sectional areax
3 (m2)
1
×
1
0
-
4
Upper limit of stress in nondesignable domainσ
¯ (MPa)
80
Upper limit of design domain’s volumeV
¯ (m3)
2.5
×
1
0
-
4
Lower limit of design variablesx
_ (m2)
5
×
1
0
-
5
Upper limit of design variablesx
¯ (m2)
2
×
1
0
-
4At the first place, the mechanical load is fixed withf
m
x
=
f
m
y
=
7000
N. Structural compliance with different temperature increases is plotted by sets of contours in Figure 3. Thick solid line and dashed line denote stress and volume constraints, respectively. Along with the lower and upper limits of the design variables, these constraints define the feasible regions of optimization, which is indicated as the dark area in the figure. The black spot denotes the optimal point where the compliance reaches a minimum. Those symbol signs are applied to Figures 3 and 4.Figure 3
Optimizations under different temperature increases.
(a)
Δ
T
=
16
°
C, f
m
x
=
f
m
y
=
7000
N
(b)
Δ
T
=
18
°
C, f
m
x
=
f
m
y
=
7000
N
(c)
Δ
T
=
20
°
C, f
m
x
=
f
m
y
=
7000
N
(d)
Δ
T
=
25
°
C, f
m
x
=
f
m
y
=
7000
NFigure 4
Optimizations under different mechanical loads.
(a)
Δ
T
=
20
°
C, f
m
x
=
f
m
y
=
5500
N
(b)
Δ
T
=
20
°
C, f
m
x
=
f
m
y
=
6500
N
(c)
Δ
T
=
20
°
C, f
m
x
=
f
m
y
=
7500
N
(d)
Δ
T
=
20
°
C, f
m
x
=
f
m
y
=
7800
NAs illustrated in Figure3, when the structure undergoes thermomechanical loads, a small temperature increase will lead to nonmonotone compliance with respect to the design variables. Consequently, the volume constraint may not be active as shown in Figures 3(b) to 3(d). In this case, using more material could weaken the structural mechanical performance due to the thermomechanical loading condition.In another case, fixingΔ
T
=
20
°C, the mechanical loads are proportionally raised. The contour lines of compliance are shown in Figure 4. At the beginning, the stress is primarily caused by thermal load. Very compliant structures are used to offset thermal stress. Typically, when both mechanical and thermal loads are small as in Figure 4(a), the optimal point lies at the lower left corner of the feasible region. At this point, the stress constraint is inactive. As mechanical loads increase, more materials are required to strengthen the structure. The volume constraint is finally active in Figure 4. We can foresee that as the mechanical loads increase, there will be no feasible design for this problem.Notably, there exist some critical optimal points for thermomechanical loads, for example, in Figures3(a) and 4(d). The volume constraint is exactly active there. A slight increase of temperature or decrease of mechanical loads will lead to inactive volume constraint. These phenomena also happen in topology optimization of continuum structures with thermomechanical loads, where inactive volume constraints and elements with intermediate densities are always found (see [20]).


## 3. Formulation of Thermomechanical Topology Optimization

Based on the analyses in the previous section, we propose to use the topology optimization formulation as(6)
Find:
X
=
x
1
,
x
2
,
…
,
x
n
Min:
C
=
1
2
F
T
U
=
1
2
F
t
h
+
F
m
T
U
s.t.:
V
X
=
∑
i
=
1
n
x
i
v
i
≤
V
¯
;
m
a
x
j
=
1,2
,
…
,
m
σ
j
V
M
≤
σ
¯
;
0
<
δ
≤
x
i
≤
1
,
i
=
1,2
,
…
,
n
,where n and m are the numbers of elements in elastic support and nondesignable domain, respectively. X is the vector of design variables, that is, the pseudodensities which vary between 0 and 1 to describe the material distribution over the design domain. C is the global compliance of the structure. F
t
h and F
m are the thermal and mechanical load vectors composing the nodal load vector F. σ
j
V
M is the Von-Mises stress of the j
t
h element in the nondesignable domain. V
(
X
) is the volume of the design domain and v
i is the volume of the i
t
h element. δ is a small constant set as 0.001 in this paper to avoid singularity of stiffness matrix.


### 3.1. Design-Dependent Thermal Load

Unlike the static mechanical load, thermal load is a typical design-dependent load depending upon the material layout over the design domain. Thus, the matching relation between the thermal load and stiffness should be carefully handled.According to the thermoelastic theory, the thermal load over theith element is expressed as(7)
F
i
t
h
=
∫
Ω
i
B
i
T
D
i
ε
i
t
h
d
Ω
,where B
i and D
i are the element strain-displacement matrix and elasticity matrix, respectively. ε
i
t
h is the thermal strain of the ith element. D
i is evaluated in terms of Young’s modulus E
i. Here, a polynomial interpolation is used to link the pseudodensity variables to element elastic modulus [23]:(8)
E
i
=
1
-
β
x
i
p
+
β
x
i
E
i
0
,where p
=
3 and β
=
1
/
16. This polynomial model works well in topology optimization with design-dependent loads, such as inertial loads, dynamic loads, and thermal loads. It maintains a positive gradient when the pseudodensity is zero.Suppose the thermal expansion coefficient is temperature-independent and only the steady-state temperature field is taken into account. The thermal strain vector can be written as follows:(9)
ε
i
t
h
=
α
i
Δ
T
i
α
T
.Here, α
i is the thermal expansion coefficient. Δ
T
i denotes the temperature increase. α is a constant vector for the calculation of strain vector which is [
1,1
,
1,0
,
0,0
] in this paper [17].Referring to the conception of thermal stress coefficient (TSC) [14], the element stiffness and the thermal stress load should be penalized independently in terms of element pseudodensity. A linear interpolation to the thermal expansion coefficient reads(10)
α
i
=
x
i
α
i
0
,where α
i
0 is the original thermal expansion coefficient. Consequently, the elemental thermal expansion coefficient now varies with the pseudodensities in coordination with Young’s modulus.With the substitution of (8)–(10) into (7), the thermal load can be expressed as(11)
F
i
t
h
=
γ
i
x
i
F
i
0
t
h
,
F
i
0
t
h
=
∫
Ω
i
B
i
T
D
i
0
α
d
Ω
,
γ
i
=
E
i
α
i
,where F
i
0
t
h is the thermal load when the element is solid with unit Young’s modulus, unit thermal expansion coefficient, and unit temperature increase which is always constant in the optimization process. γ
i denotes the thermal stress coefficient (TSC) which can be treated as an inherent material property. With the introduction of TSC, the thermal load can be explicitly expressed as the function of design variables.


### 3.2. The Global Stress Measure

Difficulties of min-max problems are generally involved in stress-based topology optimization. Naturally, stress is evaluated and constrained in each element. In this way, the overall stress level is restricted below the prescribed limit. But the scale of the optimization problem and the computational cost increase dramatically. To remedy the drawback, aggregation methods combining all local stresses into one or several Kreisselmeier-Steinhauser (KS) orp-norm global functions [24, 25] are generally used. A tradeoff between the global and local stress approaches is that global stress measure reduces the computational cost while it introduces high nonlinearity in global measure functions and loses the control over local stresses. Different aggregation strategies have been brought in. For example, the active set approach [26] aggregates the active constraints into one set, which highly reduces the number of stress constraints. However, the changes of the active set can lead to poor convergence. Le et al. [27] introduced an adaptively updated p-norm aggregation method, which can precisely approximate the maximum local stress. A clustered method [28] was proposed as a compromise between global and local methods. The constraints were aggregated into several clusters. This method was proven to be efficient and robust by numerical examples but the results strongly depended upon the number of clusters.Concerning the optimization of elastic support structure in this paper, the thermal stress constraints are imposed only in the nondesignable domain. Compared with the total elements number, the number of evaluated stresses is comparatively small. Without loss of generality, the standardp-norm measure is utilized which reads(12)
σ
P
N
=
∑
j
=
1
m
σ
j
V
M
P
1
/
P
,where σ
P
N is the p-norm of the elemental Von-Mises stress and P is the aggregation parameter. Apparently, the p-norm degenerates to the summation of the stress components when P is assigned as 1. Moreover, it can be proven mathematically that(13)
lim
p
→
∞
⁡
∑
j
=
1
m
σ
j
V
M
P
1
/
P
=
max
j
=
1,2
,
…
,
m
⁡
σ
j
V
M
.Obviously, with a sufficiently large value for P, the global measure can exactly match the maximum stress. However, the large value may lead to oscillation and ill-posed problems [29]. Consequently, to balance the effect of constraints’ quality and convergence efficiency, the value of P should be properly chosen, which is assigned as 15 in this paper.


### 3.3. Adjoint Sensitivity Analysis

Firstly, the sensitivity of the global compliance with the design-dependent thermal load is derived and written as(14)
∂
C
∂
x
i
=
1
2
U
T
∂
F
t
h
∂
x
i
+
1
2
F
m
+
F
t
h
T
·
∂
U
∂
x
i
=
1
2
U
T
∂
F
t
h
∂
x
i
+
1
2
F
m
+
F
t
h
T
·
K
-
1
∂
F
t
h
∂
x
i
-
∂
K
∂
x
i
U
=
U
T
∂
F
t
h
∂
x
i
-
1
2
U
T
∂
K
∂
x
i
U
.Sensitivities of the thermal load with respect to the design variables can be obtained by differentiating its definition formulated in (11). The derivative of the stiffness matrix can be easily evaluated according to the interpolation functions in (8).Sensitivity analysis is then carried out to evaluate the global stress measure in the nondesignable domain with respect to the pseudodensity variables, which can be written as(15)
∂
σ
P
N
∂
x
i
=
∑
j
=
1
m
σ
j
V
M
P
1
/
P
-
1
∑
j
=
1
m
σ
j
V
M
P
-
1
·
∂
σ
j
V
M
∂
x
i
=
σ
P
N
1
-
P
∑
j
=
1
m
σ
j
V
M
P
-
1
·
∂
σ
j
V
M
∂
x
i
.Here, the Von-Mises equivalent stress is defined as(16)
σ
j
V
M
σ
j
T
V
σ
j
=
B
j
U
j
-
ε
j
t
h
T
D
j
T
V
D
j
B
j
U
j
-
ε
j
t
h
,where σ
j is the element stress vector. U
j and ε
j
t
h are the element nodal displacement vector and thermal strain vector. The matrix V is a constant matrix; that is,(17)
V
2
D
=
1
-
1
2
0
-
1
2
1
0
0
0
3
,
V
3
D
=
2
-
1
-
1
-
1
2
-
1
-
1
-
1
2
0
3
×
3
0
3
×
3
6
6
6
.As only the Von-Mises stresses in the nondesignable domain are concerned, we have(18)
∂
σ
j
V
M
∂
x
i
=
B
j
U
j
-
ε
j
t
h
T
D
j
T
V
D
j
B
j
U
j
-
ε
j
t
h
T
D
j
T
V
D
j
B
j
U
j
-
ε
j
t
h
B
j
∂
U
j
∂
x
i
=
1
σ
j
V
M
·
σ
j
T
V
D
j
B
j
∂
U
j
∂
x
i
,where ∂
D
j
/
∂
x
i
=
∂
ε
j
t
h
/
∂
x
i
=
0. In this case, the singularity phenomenon in stress-based topology optimization is actually not involved in this optimization problem.We can also define(19)
U
j
=
A
j
U
,
∂
U
j
∂
x
i
=
A
j
∂
U
∂
x
i
.The matrix A transforms the global nodal displacement vector to an element one.Based on the differentiation of finite element equilibrium equation, we have the derivative of the nodal displacement vector(20)
∂
U
∂
x
i
=
K
-
1
∂
F
t
h
∂
x
i
-
∂
K
∂
x
i
U
.The substitution of the above equation into the derivative of p-norm aggregation function produces(21)
∂
σ
P
N
∂
x
i
=
σ
P
N
1
-
P
∑
j
=
1
m
σ
j
V
M
P
-
2
σ
j
T
V
D
B
j
A
j
K
-
1
∂
F
t
h
∂
x
i
-
∂
K
∂
x
i
U
.Suppose(22)
Q
T
=
σ
P
N
1
-
P
∑
j
=
1
m
σ
j
V
M
P
-
2
σ
j
T
V
D
B
j
A
j
K
-
1is adjoint displacement; the derivative of the p-norm function can be obtained by an additional finite element calculation.


## 3.1. Design-Dependent Thermal Load

Unlike the static mechanical load, thermal load is a typical design-dependent load depending upon the material layout over the design domain. Thus, the matching relation between the thermal load and stiffness should be carefully handled.According to the thermoelastic theory, the thermal load over theith element is expressed as(7)
F
i
t
h
=
∫
Ω
i
B
i
T
D
i
ε
i
t
h
d
Ω
,where B
i and D
i are the element strain-displacement matrix and elasticity matrix, respectively. ε
i
t
h is the thermal strain of the ith element. D
i is evaluated in terms of Young’s modulus E
i. Here, a polynomial interpolation is used to link the pseudodensity variables to element elastic modulus [23]:(8)
E
i
=
1
-
β
x
i
p
+
β
x
i
E
i
0
,where p
=
3 and β
=
1
/
16. This polynomial model works well in topology optimization with design-dependent loads, such as inertial loads, dynamic loads, and thermal loads. It maintains a positive gradient when the pseudodensity is zero.Suppose the thermal expansion coefficient is temperature-independent and only the steady-state temperature field is taken into account. The thermal strain vector can be written as follows:(9)
ε
i
t
h
=
α
i
Δ
T
i
α
T
.Here, α
i is the thermal expansion coefficient. Δ
T
i denotes the temperature increase. α is a constant vector for the calculation of strain vector which is [
1,1
,
1,0
,
0,0
] in this paper [17].Referring to the conception of thermal stress coefficient (TSC) [14], the element stiffness and the thermal stress load should be penalized independently in terms of element pseudodensity. A linear interpolation to the thermal expansion coefficient reads(10)
α
i
=
x
i
α
i
0
,where α
i
0 is the original thermal expansion coefficient. Consequently, the elemental thermal expansion coefficient now varies with the pseudodensities in coordination with Young’s modulus.With the substitution of (8)–(10) into (7), the thermal load can be expressed as(11)
F
i
t
h
=
γ
i
x
i
F
i
0
t
h
,
F
i
0
t
h
=
∫
Ω
i
B
i
T
D
i
0
α
d
Ω
,
γ
i
=
E
i
α
i
,where F
i
0
t
h is the thermal load when the element is solid with unit Young’s modulus, unit thermal expansion coefficient, and unit temperature increase which is always constant in the optimization process. γ
i denotes the thermal stress coefficient (TSC) which can be treated as an inherent material property. With the introduction of TSC, the thermal load can be explicitly expressed as the function of design variables.


## 3.2. The Global Stress Measure

Difficulties of min-max problems are generally involved in stress-based topology optimization. Naturally, stress is evaluated and constrained in each element. In this way, the overall stress level is restricted below the prescribed limit. But the scale of the optimization problem and the computational cost increase dramatically. To remedy the drawback, aggregation methods combining all local stresses into one or several Kreisselmeier-Steinhauser (KS) orp-norm global functions [24, 25] are generally used. A tradeoff between the global and local stress approaches is that global stress measure reduces the computational cost while it introduces high nonlinearity in global measure functions and loses the control over local stresses. Different aggregation strategies have been brought in. For example, the active set approach [26] aggregates the active constraints into one set, which highly reduces the number of stress constraints. However, the changes of the active set can lead to poor convergence. Le et al. [27] introduced an adaptively updated p-norm aggregation method, which can precisely approximate the maximum local stress. A clustered method [28] was proposed as a compromise between global and local methods. The constraints were aggregated into several clusters. This method was proven to be efficient and robust by numerical examples but the results strongly depended upon the number of clusters.Concerning the optimization of elastic support structure in this paper, the thermal stress constraints are imposed only in the nondesignable domain. Compared with the total elements number, the number of evaluated stresses is comparatively small. Without loss of generality, the standardp-norm measure is utilized which reads(12)
σ
P
N
=
∑
j
=
1
m
σ
j
V
M
P
1
/
P
,where σ
P
N is the p-norm of the elemental Von-Mises stress and P is the aggregation parameter. Apparently, the p-norm degenerates to the summation of the stress components when P is assigned as 1. Moreover, it can be proven mathematically that(13)
lim
p
→
∞
⁡
∑
j
=
1
m
σ
j
V
M
P
1
/
P
=
max
j
=
1,2
,
…
,
m
⁡
σ
j
V
M
.Obviously, with a sufficiently large value for P, the global measure can exactly match the maximum stress. However, the large value may lead to oscillation and ill-posed problems [29]. Consequently, to balance the effect of constraints’ quality and convergence efficiency, the value of P should be properly chosen, which is assigned as 15 in this paper.


## 3.3. Adjoint Sensitivity Analysis

Firstly, the sensitivity of the global compliance with the design-dependent thermal load is derived and written as(14)
∂
C
∂
x
i
=
1
2
U
T
∂
F
t
h
∂
x
i
+
1
2
F
m
+
F
t
h
T
·
∂
U
∂
x
i
=
1
2
U
T
∂
F
t
h
∂
x
i
+
1
2
F
m
+
F
t
h
T
·
K
-
1
∂
F
t
h
∂
x
i
-
∂
K
∂
x
i
U
=
U
T
∂
F
t
h
∂
x
i
-
1
2
U
T
∂
K
∂
x
i
U
.Sensitivities of the thermal load with respect to the design variables can be obtained by differentiating its definition formulated in (11). The derivative of the stiffness matrix can be easily evaluated according to the interpolation functions in (8).Sensitivity analysis is then carried out to evaluate the global stress measure in the nondesignable domain with respect to the pseudodensity variables, which can be written as(15)
∂
σ
P
N
∂
x
i
=
∑
j
=
1
m
σ
j
V
M
P
1
/
P
-
1
∑
j
=
1
m
σ
j
V
M
P
-
1
·
∂
σ
j
V
M
∂
x
i
=
σ
P
N
1
-
P
∑
j
=
1
m
σ
j
V
M
P
-
1
·
∂
σ
j
V
M
∂
x
i
.Here, the Von-Mises equivalent stress is defined as(16)
σ
j
V
M
σ
j
T
V
σ
j
=
B
j
U
j
-
ε
j
t
h
T
D
j
T
V
D
j
B
j
U
j
-
ε
j
t
h
,where σ
j is the element stress vector. U
j and ε
j
t
h are the element nodal displacement vector and thermal strain vector. The matrix V is a constant matrix; that is,(17)
V
2
D
=
1
-
1
2
0
-
1
2
1
0
0
0
3
,
V
3
D
=
2
-
1
-
1
-
1
2
-
1
-
1
-
1
2
0
3
×
3
0
3
×
3
6
6
6
.As only the Von-Mises stresses in the nondesignable domain are concerned, we have(18)
∂
σ
j
V
M
∂
x
i
=
B
j
U
j
-
ε
j
t
h
T
D
j
T
V
D
j
B
j
U
j
-
ε
j
t
h
T
D
j
T
V
D
j
B
j
U
j
-
ε
j
t
h
B
j
∂
U
j
∂
x
i
=
1
σ
j
V
M
·
σ
j
T
V
D
j
B
j
∂
U
j
∂
x
i
,where ∂
D
j
/
∂
x
i
=
∂
ε
j
t
h
/
∂
x
i
=
0. In this case, the singularity phenomenon in stress-based topology optimization is actually not involved in this optimization problem.We can also define(19)
U
j
=
A
j
U
,
∂
U
j
∂
x
i
=
A
j
∂
U
∂
x
i
.The matrix A transforms the global nodal displacement vector to an element one.Based on the differentiation of finite element equilibrium equation, we have the derivative of the nodal displacement vector(20)
∂
U
∂
x
i
=
K
-
1
∂
F
t
h
∂
x
i
-
∂
K
∂
x
i
U
.The substitution of the above equation into the derivative of p-norm aggregation function produces(21)
∂
σ
P
N
∂
x
i
=
σ
P
N
1
-
P
∑
j
=
1
m
σ
j
V
M
P
-
2
σ
j
T
V
D
B
j
A
j
K
-
1
∂
F
t
h
∂
x
i
-
∂
K
∂
x
i
U
.Suppose(22)
Q
T
=
σ
P
N
1
-
P
∑
j
=
1
m
σ
j
V
M
P
-
2
σ
j
T
V
D
B
j
A
j
K
-
1is adjoint displacement; the derivative of the p-norm function can be obtained by an additional finite element calculation.


## 4. Numerical Study

In this section, a 2D cantilever beam and a 3D nozzle flap of a turbine engine are tested to study the topology optimization of elastic supports under thermomechanical loads.


### 4.1. 2D Cantilever Beam

As shown in Figure5, a 2D cantilever beam consists of a nondesignable top surface and an elastic support assigned as design domain. The material properties are listed in Table 2. Since the thermal expansion coefficient of elastic support is higher than that of the nondesignable domain, significant thermal stresses will occur as a result of temperature increase.Table 2
Material properties used in the optimization.
Properties
Elastic support
Nondesignable domain
Elastic modulus (MPa)
71000
3500
Coefficient of thermal expansion (10−6/°C)
23
7
Poisson’s ratio
0.33
0.4Figure 5
Definition of the 2D cantilever beam.In accordance with Section3, the global compliance of the whole structure is minimized as the design objective. The volume fraction is constrained to 50% and the upper bound of the stress constraints in the nondesignable domain is 200 MPa.The optimized results with various load combinations are listed and compared in Figure6. The result in Figure 6(a) is set as a benchmark with clear configuration. Firstly, the mechanical load is fixed at 20 N/mm in Figures 6(a)–6(d). As the temperature increases, the volume constraint starts to be inactive. More elements with intermediate densities arise. This is due to the nature of thermal stress constraint where the elements with compliant materials can offset the thermal stress better. Later, considering the changes in Figures 6(c), 6(e), and 6(f), clear structure patterns gradually appear when larger mechanical loads are applied. Thus, the elastic support is able to undergo the pressure with good stiffness. However, when the thermomechanical loads are too high, the optimization finds no feasible solution where the stress constraints are violated as shown in Figure 6(f).Figure 6
Optimized results of minimizing compliance under different loads.
(a)
Δ
T
=
50
°C, q = 20N/mm. Optimized volume fraction 50%. Maximum Stress 154 MPa
(b)
Δ
T
=
200
°C, q = 20N/mm. Optimized volume fraction 50%. Maximum stress 166.4 MPa
(c)
Δ
T
=
500
°C, q = 20N/mm. Optimized volume fraction 41.7%. Maximum stress 200 MPa
(d)
Δ
T
=
1000
°C, q = 20N/mm. Optimized volume fraction 35.8%. Maximum stress 200 MPa
(e)
Δ
T
=
500
°C, q = 40N/mm. Optimized volume fraction 50%. Maximum stress 200 MPa
(f)
Δ
T
=
500
°C, q = 200N/mm. Optimized volume fraction 50%. Maximum stress 812.4 MPaThe optimized results have shown that the relative magnitude of the thermal and mechanical load greatly influences the optimized results, which have good consistency with the analytical solution of the three-bar truss model in Section3. Moreover, the existence of a large amount of intermediate density material is reasonable when the thermal load is dominant over the mechanical load. Eliminating the grey elements directly or using some numerical schemes may improve the global stiffness but will unfortunately lead to higher stress level.


### 4.2. 3D Nozzle Flap

As shown in Figure7, a nozzle flap of a turbine engine is composed of titanium stiffeners and a ceramic plate. Material properties are listed in Table 3. Figure 7 also illustrates the thermoelastic loads applied on the model, including a uniform pressure of 1 MPa and global temperature increase of 500°C. Significant thermal stresses are generated in the plate due to the different thermal expansion coefficient. In practical design, the stiffeners are assigned as the elastic support design domain and the plate is nondesignable. The design objective is to minimize the global compliance with a 25 MPa stress constraint on the nondesignable plate and a 30% material volume constraint on the design domain.Table 3
Material properties used in the optimization.
Material properties
Titanium
Ceramic
Elastic modulus (MPa)
11000
3500
Coefficient of thermal expansion (10−6/°C)
10
7
Poisson’s ratio
0.33
0.4Figure 7
Nozzle flap model and its design domain, loads, and boundary conditions.The optimized design of elastic support as shown in Figure8 is obtained by topology optimization using the proposed formulation. The optimized design is presented by hiding the elements with their pseudodensities under 0.5 to show a clear structural configuration. The two strongest stiffeners are composed of solid elements while the detailed structural branches are using intermediate material with pseudodensities between 0.5 and 0.9.Figure 8
Topological optimized design with a temperature increase of 500°C.To further emphasize the effect of the thermal loads, two more designs are obtained by using different temperature increases, that is, 100°C and 1000°C, with identical mechanical loads. The optimized designs are shown in Figure9. Compared with the structural topology in Figure 8, the optimized design in Figure 9(a) has shown a much clearer load carrying path with less intermediate material as the mechanical load is dominant. In Figure 9(b), an extremely high temperature increase of 1000°C is used. The optimized design is mostly composed of intermediate material as expected. No clear structural configuration is achieved.Figure 9
Topological optimized design with different temperature increases.
(a)
100°C
(b)
1000°CCAD model as shown in Figure10(a) is then rebuilt according to the optimized design in Figure 8. An existing design of the nozzle flap is shown in Figure 10(b) for comparison. To verify the effect of topology optimization, the two models are analyzed, respectively, with refined finite element mesh. Two models share identical boundary conditions and thermomechanical loads; the stress distribution in the bottom plate is shown in Figure 11. The overall comparison of the two designs is shown in Table 4.Table 4
Comparison of the optimal design and the existing design.
Topological optimized design
Rebuilt optimized design
Existing design
Maximum stress in the bottom plate (MPa)
25.00
25.00
284.7
Compliance of the whole structure (KJ)
302.2
254.6
271.9
Volume of the elastic support (mm3)
2.482 × 107
2.078 × 107
2.207 × 107Figure 10
Rebuilt optimized design and existing design to be compared.
(a)
Rebuilt optimized design
(b)
Existing designFigure 11
Stress distribution in the bottom plate (MPa).
(a)
Optimized design
(b)
Existing designCompared with the existing design, the optimized design reduces the maximum stress in the bottom plate significantly from 284.7 MPa to 25 MPa. The global compliance decreases from 271.9 KJ to 254.6 KJ. Meanwhile, material of 0.129 × 107 mm3 is saved, that is, 5.676 kg lighter than before.


## 4.1. 2D Cantilever Beam

As shown in Figure5, a 2D cantilever beam consists of a nondesignable top surface and an elastic support assigned as design domain. The material properties are listed in Table 2. Since the thermal expansion coefficient of elastic support is higher than that of the nondesignable domain, significant thermal stresses will occur as a result of temperature increase.Table 2
Material properties used in the optimization.
Properties
Elastic support
Nondesignable domain
Elastic modulus (MPa)
71000
3500
Coefficient of thermal expansion (10−6/°C)
23
7
Poisson’s ratio
0.33
0.4Figure 5
Definition of the 2D cantilever beam.In accordance with Section3, the global compliance of the whole structure is minimized as the design objective. The volume fraction is constrained to 50% and the upper bound of the stress constraints in the nondesignable domain is 200 MPa.The optimized results with various load combinations are listed and compared in Figure6. The result in Figure 6(a) is set as a benchmark with clear configuration. Firstly, the mechanical load is fixed at 20 N/mm in Figures 6(a)–6(d). As the temperature increases, the volume constraint starts to be inactive. More elements with intermediate densities arise. This is due to the nature of thermal stress constraint where the elements with compliant materials can offset the thermal stress better. Later, considering the changes in Figures 6(c), 6(e), and 6(f), clear structure patterns gradually appear when larger mechanical loads are applied. Thus, the elastic support is able to undergo the pressure with good stiffness. However, when the thermomechanical loads are too high, the optimization finds no feasible solution where the stress constraints are violated as shown in Figure 6(f).Figure 6
Optimized results of minimizing compliance under different loads.
(a)
Δ
T
=
50
°C, q = 20N/mm. Optimized volume fraction 50%. Maximum Stress 154 MPa
(b)
Δ
T
=
200
°C, q = 20N/mm. Optimized volume fraction 50%. Maximum stress 166.4 MPa
(c)
Δ
T
=
500
°C, q = 20N/mm. Optimized volume fraction 41.7%. Maximum stress 200 MPa
(d)
Δ
T
=
1000
°C, q = 20N/mm. Optimized volume fraction 35.8%. Maximum stress 200 MPa
(e)
Δ
T
=
500
°C, q = 40N/mm. Optimized volume fraction 50%. Maximum stress 200 MPa
(f)
Δ
T
=
500
°C, q = 200N/mm. Optimized volume fraction 50%. Maximum stress 812.4 MPaThe optimized results have shown that the relative magnitude of the thermal and mechanical load greatly influences the optimized results, which have good consistency with the analytical solution of the three-bar truss model in Section3. Moreover, the existence of a large amount of intermediate density material is reasonable when the thermal load is dominant over the mechanical load. Eliminating the grey elements directly or using some numerical schemes may improve the global stiffness but will unfortunately lead to higher stress level.


## 4.2. 3D Nozzle Flap

As shown in Figure7, a nozzle flap of a turbine engine is composed of titanium stiffeners and a ceramic plate. Material properties are listed in Table 3. Figure 7 also illustrates the thermoelastic loads applied on the model, including a uniform pressure of 1 MPa and global temperature increase of 500°C. Significant thermal stresses are generated in the plate due to the different thermal expansion coefficient. In practical design, the stiffeners are assigned as the elastic support design domain and the plate is nondesignable. The design objective is to minimize the global compliance with a 25 MPa stress constraint on the nondesignable plate and a 30% material volume constraint on the design domain.Table 3
Material properties used in the optimization.
Material properties
Titanium
Ceramic
Elastic modulus (MPa)
11000
3500
Coefficient of thermal expansion (10−6/°C)
10
7
Poisson’s ratio
0.33
0.4Figure 7
Nozzle flap model and its design domain, loads, and boundary conditions.The optimized design of elastic support as shown in Figure8 is obtained by topology optimization using the proposed formulation. The optimized design is presented by hiding the elements with their pseudodensities under 0.5 to show a clear structural configuration. The two strongest stiffeners are composed of solid elements while the detailed structural branches are using intermediate material with pseudodensities between 0.5 and 0.9.Figure 8
Topological optimized design with a temperature increase of 500°C.To further emphasize the effect of the thermal loads, two more designs are obtained by using different temperature increases, that is, 100°C and 1000°C, with identical mechanical loads. The optimized designs are shown in Figure9. Compared with the structural topology in Figure 8, the optimized design in Figure 9(a) has shown a much clearer load carrying path with less intermediate material as the mechanical load is dominant. In Figure 9(b), an extremely high temperature increase of 1000°C is used. The optimized design is mostly composed of intermediate material as expected. No clear structural configuration is achieved.Figure 9
Topological optimized design with different temperature increases.
(a)
100°C
(b)
1000°CCAD model as shown in Figure10(a) is then rebuilt according to the optimized design in Figure 8. An existing design of the nozzle flap is shown in Figure 10(b) for comparison. To verify the effect of topology optimization, the two models are analyzed, respectively, with refined finite element mesh. Two models share identical boundary conditions and thermomechanical loads; the stress distribution in the bottom plate is shown in Figure 11. The overall comparison of the two designs is shown in Table 4.Table 4
Comparison of the optimal design and the existing design.
Topological optimized design
Rebuilt optimized design
Existing design
Maximum stress in the bottom plate (MPa)
25.00
25.00
284.7
Compliance of the whole structure (KJ)
302.2
254.6
271.9
Volume of the elastic support (mm3)
2.482 × 107
2.078 × 107
2.207 × 107Figure 10
Rebuilt optimized design and existing design to be compared.
(a)
Rebuilt optimized design
(b)
Existing designFigure 11
Stress distribution in the bottom plate (MPa).
(a)
Optimized design
(b)
Existing designCompared with the existing design, the optimized design reduces the maximum stress in the bottom plate significantly from 284.7 MPa to 25 MPa. The global compliance decreases from 271.9 KJ to 254.6 KJ. Meanwhile, material of 0.129 × 107 mm3 is saved, that is, 5.676 kg lighter than before.


## 5. Conclusion

In this paper, topology optimization of elastic supporting structures under thermomechanical loads is investigated. A three-bar truss model is firstly employed to reveal the particularity of thermoelastic problems, that is, nonmonotonous compliance, inactive volume constraint with high temperature increase, and so forth. Similar appearances also have been found in the topology optimization of a 2D cantilever beam structure presented in this paper. On account of stress-based topology optimization with large numbers of design constraints, global stress measure approach based onp-norm function is used to aggregate the stress constraints involved in each iteration. Meanwhile, formulation of design-dependent thermal load is presented with different interpolation functions assigned for material elastic modulus and thermal expansion coefficients. Sensitivity analysis is then carried out to evaluate the global stress measure in nondesignable domain and the compliance with respect to the pseudodensity variables. In the Numerical Study, a 2D cantilever beam model and a 3D nozzle flap are optimized. Compared with the existing design, the optimized designs not only use fewer materials, but also are both stiffer and better in reducing the maximum stress in nondesignable domain.



---

*Source: 7372603-2016-06-27.xml*