# An Enhanced Analytical Target Cascading and Kriging Model Combined Approach for Multidisciplinary Design Optimization

Authors: Ping Jiang; Jianzhuang Wang; Qi Zhou; Xiaolin Zhang
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

Multidisciplinary design optimization (MDO) has been applied widely in the design of complex engineering systems. To ease MDO problems, analytical target cascading (ATC) organizes MDO process into multilevels according to the components of engineering systems, which provides a promising way to deal with MDO problems. ATC adopts a coordination strategy to coordinate the couplings between two adjacent levels in the design optimization process; however, existing coordination strategies in ATC face the obstacles of complicated coordination process and heavy computation cost. In order to conquer this problem, a quadratic exterior penalty function (QEPF) based ATC (QEPF-ATC) approach is proposed, where QEPF is adopted as the coordination strategy. Moreover, approximate models are adopted widely to replace the expensive simulation models in MDO; a QEPF-ATC and Kriging model combined approach is further proposed to deal with MDO problems, owing to the comprehensive performance, high approximation accuracy, and robustness of Kriging model. Finally, the geometric programming and reducer design cases are given to validate the applicability and efficiency of the proposed approach.


---

## Body

## 1. Introduction

Multidisciplinary design optimization (MDO) is a methodology to design complex engineering systems by exploiting the interactions between disciplines and has been successfully applied in the design of complex engineering systems, such as aircraft, automobiles, and mechanical equipment [1, 2]. MDO method, also known as MDO strategy, is a hot topic in the research area of MDO currently [3]. Depending on whether the process of design optimization is hierarchically designed, the MDO method is divided into two categories: single-level MDO method and multilevels MDO method. Compared with single-level MDO method, the multilevels MDO method can better match the organization structure of complex engineering systems and is beneficial in implementing concurrent design and distributed computation. There are four main multilevels methods: collaborative optimization (CO) [4, 5], concurrent subspace optimization (CSSO) [6, 7], bilevel integrated system synthesis (BLISS) [8, 9], and analytical target cascading (ATC) [10, 11]. CO, CSSO, and BLISS are all two-level MDO methods. Compared with CO, CSSO, and BLISS methods, ATC can organize MDO process into different levels according to the components of engineering systems, and the convergence of ATC method has been proofed strictly [10], so ATC provides a promising way to deal with MDO problems.Currently, ATC needs to adopt a coordination strategy to coordinate the couplings between two adjacent levels in the design optimization process, and an appropriate coordination strategy can improve the computation efficiency of ATC. Till now, some coordination strategies have been adopted in ATC, such as augmented Lagrangian function (ALF) [12], diagonal quadratic approximation (DQA) [13], and Lagrangian dual function (LDF) [14]. However, these coordination strategies based ATC approaches face the obstacles of complicated coordination process and heavy computation cost in dealing with MDO problems. In order to conquer this problem, a quadratic exterior penalty function (QEPF) based ATC (QEPF-ATC) approach is proposed in this paper, where QEPF is adopted as the coordination strategy. Compared with exiting ATC method, the QEPF-ATC approach can simplify the coordination process and reduce the computation cost.In the design optimization of complex engineering systems, generally, a large number of simulation models with high accuracy (such as the structural finite element analysis models, aerodynamic analysis models, and computational fluid dynamic models) are involved. These simulation models increase the computation expense greatly. Therefore, approximate models are adopted widely to replace the expensive simulation models in MDO [15–18]. In this paper, QEPF-ATC and Kriging model [15, 18] combined approach is further proposed to deal with MDO problems, owing to the comprehensive performance, high approximation accuracy, and robustness of Kriging model. Finally, the geometric programming and reducer design cases are given to validate usability and effectiveness of the proposed QEPF-ATC approach and the QEPF-ATC and Kriging model combined approach, respectively.The rest of this paper is organized as follows. In Section2, ATC is introduced, and the QEPF-ATC approach is proposed. QEPF-ATC and Kriging model combined approach is implemented in Section 3. In Section 4, the geometric programming and reducer design cases are given to validate the usability and effectiveness of the proposed approach. Conclusion and future work are discussed in Section 5.


## 2. The Proposed QEPF-ATC Approach

### 2.1. ATC

The design optimization problem of complex engineering systems is decomposed into a set of subproblems in ATC approach, where the desired design specifications or targets are established at the top level and cascaded down to each lower level subproblem. The target transformation is a decomposition process from up to down in a level-by-level way. Moreover, subproblems need to get the responses to match the cascaded targets. If targets cascaded to subsystems are not achievable or compatible, the feedback from subsystems is required, resulting in an iterative target cascading process. Hence, targets and responses are updated and coordinated iteratively to achieve global consistency (see Figure1) [10]. Coordination strategy is needed to implement the communication of target and response values between parents and children problems.Figure 1
Schematic of the ATC process.ATC is a model-based optimization method. There exist two kinds of models, optimal design modelsP and analysis models R. The former is just so-called optimization design problem that includes optimization object, design variables, and constraints. The latter is physical model or mathematical model that is the foundation of the optimization design. Figure 2 shows the data flow between subproblems (parent and children) [19].Figure 2
Information flow between children subproblem and parent subproblem.In Figure2, there exist consistency constraints for each subproblem. In order to make each subproblem independent, consistency constraint is introduced as a kind of equality constraints. For the subproblem at a given level, consistency constraints are as follows:(1)
c
1
′
=
t
i
j
-
r
i
j
=
0
,
c
1
′
=
R
i
j
i
-
1
-
R
i
j
i
=
0
,
c
1
′
=
(
y
i
p
i
-
1
-
y
i
j
i
)
=
0
,
c
2
′
=
t
i
+
1
k
-
r
i
+
1
k
=
0
,
c
2
′
=
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
=
0
,
c
2
′
=
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
=
0
.The general mathematical notation for subproblemsP
i
j in the hierarchy is then(2)
P
i
j
:
min
⁡
f
i
j
w
.
r
.
t
.
k
x
i
j
,
y
i
j
i
,
R
(
i
+
1
)
j
i
,
y
(
i
+
1
)
j
i
where
R
i
j
=
r
i
j
(
R
(
i
+
1
)
j
i
,
x
i
j
,
y
i
j
i
)
s
.
t
.
k
c
1
=
(
R
i
j
i
-
1
-
R
i
j
i
)
=
0
,
c
2
=
(
y
i
p
i
-
1
-
y
i
j
i
)
=
0
,
k
k
c
3
=
(
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
)
=
0
,
k
k
c
4
=
(
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
)
=
0
,
k
k
g
i
j
(
R
i
+
1
k
1
i
,
…
,
R
(
i
+
1
)
k
c
i
j
i
,
x
i
j
,
y
i
j
i
)
≤
0
k
k
h
i
j
(
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
)
=
0
,where P
i
j is the problem formulation of element j at level i. p is used to designate the parent of element j. k designates children of element j. x
i
j is local variables for subproblem P
i
j. r
i
j is the function that calculates responses for P
i
j. R
i
j
i is response of P
i
j. R
i
j
i
-
1 is the target for P
i
j cascaded from the (i
-
1) level parent element. R
(
i
+
1
)
k
i
+
1 is response of the (i
+
1) level children element. R
(
i
+
1
)
j
i
=
(
R
(
i
+
1
)
k
1
i
,
…
,
R
(
i
+
1
)
k
c
i
j
i
) is the target for the (i
+
1) level children element. y
i
j
i is response of linking variable of P
i
j. y
i
p
i
-
1 is the target of linking variable for P
i
j from the (i
-
1) level parent element. y
(
i
+
1
)
k
i
+
1 is response of linking variable of the (i
+
1) level children element. y
(
i
+
1
)
j
i is the target of linking variable for the (i
+
1) level children element from P
i
j. f
i
j is the local objective for P
i
j. g
i
j is function of inequality constraints for P
i
j in negative null form. h
i
j is function of equality constraints for P
i
j in null form.The coordination strategy in ATC deals with the definition of deviation between targets and responses (namely, consistency constraint) and updating weight coefficients therein. By setting appropriate norm and weight coefficients to coordinate each subproblem, the optimal design point can be obtained. The coordination strategy is the key research in ATC, and an appropriate coordination strategy can improve the computation efficiency of ATC. In recent years, some coordination strategies in ATC have been studied.


#### 2.1.1. The Coordination Strategy Based on ALF

According to the ALF, the deviation between targets and responses is presented by augmented Lagrangian function, and weight coefficients are presented by the augmented Lagrange multiplier. The coordination way based on ALF is shown as(3)
min
⁡
f
i
j
+
λ
i
j
R
(
R
i
j
i
-
R
i
j
i
-
1
)
+
w
i
j
R
∘
R
i
j
i
-
R
i
j
i
-
1
2
+
λ
i
j
y
(
y
i
p
i
-
1
-
y
i
j
i
)
+
w
i
j
y
∘
y
i
p
i
-
1
-
y
i
j
i
2
+
λ
i
+
1
j
R
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
+
w
(
i
+
1
)
j
R
∘
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
2
+
λ
i
+
1
j
y
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
+
w
i
+
1
j
y
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
2
,
s
.
t
.
g
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
≤
0
k
h
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
=
0
.Weight coefficientsλ
i
j
R, λ
i
j
y, w
i
j
R, and w
i
j
y in (3) are updated in accordance with the following equation:(4)
λ
i
j
R
k
+
1
=
λ
i
j
R
k
+
2
w
i
j
R
k
R
i
j
i
k
-
R
i
j
i
-
1
k
,
λ
i
j
y
k
+
1
=
λ
i
j
y
k
+
2
w
i
j
y
k
y
i
j
i
k
-
y
i
j
i
-
1
k
,
w
i
j
R
k
+
1
=
β
w
i
j
R
k
,
w
i
j
y
k
+
1
=
β
w
i
j
y
k
,where the superscript k indicates iterations and β is step length. Usually, β
≥
1, but 2
<
β
<
3 can accelerate convergence.


#### 2.1.2. The Coordination Strategy Based on DQA

DQA is realized by approximating the deviation between targets and responses on the basis of ALF. The specific approximation is(5)
R
i
j
i
-
R
i
j
i
-
1
2
=
R
i
j
i
2
-
2
R
i
j
i
R
i
j
i
-
1
+
R
i
j
i
-
1
2
where
,
R
i
j
i
R
i
j
i
-
1
≅
R
i
j
i
k
-
1
R
i
j
i
-
1
k
-
1
+
R
i
j
i
k
-
1
R
i
j
i
-
R
i
j
i
-
1
k
-
1
+
R
i
j
i
-
1
k
-
1
R
i
j
i
-
R
i
j
i
k
-
1
so
,
R
i
j
i
-
R
i
j
i
-
1
2
=
R
i
j
i
-
1
k
-
1
-
R
i
j
i
2
+
R
i
j
i
-
R
i
j
i
-
1
k
-
1
2
+
R
i
j
i
k
-
1
R
i
j
i
-
1
k
-
1
(
constant
value
)
.Based on DQA, the cross term among the deviation is approximated; that is to say, an order Taylor type expansion is conducted in the last design point. Among them, the weight coefficient is updated according to (5) too.


#### 2.1.3. The Coordination Strategy Based on LDF

According to the LDF, the deviation is dealt with by Lagrangian duality function, and weight coefficients are presented by the dual Lagrange multiplier. The coordination of LDF is shown as(6)
min
⁡
f
i
j
+
λ
i
j
R
R
i
j
i
-
R
i
j
i
-
1
+
μ
i
j
R
∘
R
i
j
i
-
R
i
j
i
-
1
2
kkkkkkk
+
λ
i
j
y
y
i
p
i
-
1
-
y
i
j
i
+
μ
i
j
y
∘
y
i
p
i
-
1
-
y
i
j
i
2
s
.
t
.
g
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
≤
0
k
k
k
k
l
h
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
=
0
.Weight coefficientsλ
i
j
R, λ
i
j
y, μ
i
j
R, and μ
i
j
y in (6) are updated in accordance with the following equation:(7)
λ
i
j
R
k
+
1
=
λ
i
j
R
k
+
α
i
j
R
k
R
i
j
i
k
-
R
i
j
i
-
1
k
,
λ
i
j
y
k
+
1
=
λ
i
j
y
k
+
α
i
j
y
k
y
i
j
i
k
-
y
i
j
i
-
1
k
,
μ
i
j
R
k
+
1
=
λ
i
j
R
k
,
μ
i
j
y
k
+
1
=
λ
i
j
y
k
,
α
i
j
R
k
+
1
=
1
+
m
k
+
m
·
1
R
i
j
i
k
-
R
i
j
i
-
1
k
,
α
i
j
y
k
+
1
=
1
+
m
k
+
m
·
1
y
i
j
i
k
-
y
i
j
i
-
1
k
,where the superscript k indicates iterations and m is a positive integers.The difference of these coordination strategies is the way of dealing with consistency constraints and weight coefficients. There exist multiweight coefficients in the above three coordination strategies to increase computational cost. For this reason, the coordination strategy based on quadratic exterior penalty function among which only one weight coefficient exists is put forward, which simplifies the expression of optimization problem greatly, thus improving the efficiency to make ATC have wider application. The coordination strategy based on the quadratic exterior penalty function will be described briefly.


### 2.2. The QEPF-ATC Approach

A general equality-constrained optimization problem can be defined as follows:(8)
min
⁡
f
(
x
)
s
.
t
.
lk
c
i
(
x
)
=
0
,
i
=
1,2
,
…
,
m
k
k
k
k
l
k
x
lb
≤
x
≤
x
up
,where x is the vector of the design variables with x
lb and x
up as its lower and upper bounds, respectively. The function f
(
x
) represents the objective function. c
i
(
x
) is the equality constraints.According to QEPF, consistency constraints are added into the objective function by introducing a penalty factor, and the penalty factor is used as the weighting coefficient to realize the coordination. The quadratic exterior penalty function is defined as(9)
F
x
,
ρ
=
f
x
+
ρ
2
∑
i
=
1
m
c
i
x
2
=
f
x
+
ρ
2
c
x
T
c
x
,where the weight coefficient (penalty factor) is updated in accordance with (10). The value of λ influences the convergence efficiency and computational efficiency so that it is necessary to set up appropriate λ to obtain much better solution efficiency. Consider (10)
ρ
k
+
1
=
λ
ρ
k
,
λ
=
1.4

~


10
.The mathematical formulation of QEPF-ATC is(11)
P
i
j
:
min
⁡
f
i
j
+
ρ
2
R
i
j
i
-
1
-
R
i
j
i
2
+
y
i
p
i
-
1
-
y
i
j
i
2
k
+
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
2
k
+
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
2
s
.
t
.
g
i
j
≤
0
,
h
i
j
=
0
l
R
i
j
=
r
i
j
x
i
j
¯
,
l
x
i
j
¯
=
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
.Generally, the much bigger the penalty factor is the much closer the optimization results are to the optimal solution of the original problem. But, in the practical application, the value of the penalty factor must be moderate when considering the computational cost, due to the fact that if the penalty factor is too small, the optimization results may get bigger error from the original problem to bring low computational efficiency. Otherwise, too large penalty factor would increase the computational expense.Figure3 shows the solving process of the QEPF-ATC approach, and the main steps are as follows:(1)
decompose the original problem and confirm design variables and constraints of system and subsystems;
(2)
initialize design variables (including local design variables and linking variables) andρ
(
0
) and transfer the corresponding initial value to the subsystem as the target;
(3)
solve each subsystem and pass back current results to system level;
(4)
optimize the system problem on the basis of the third step and transfer new target to each subsystem;
(5)
obtain the optimization result if the convergence condition is met; otherwise, update the weight coefficientρ according to (10) and go back to the third step until the convergence condition is met.Figure 3
The solving process of QEPF-ATC approach.


## 2.1. ATC

The design optimization problem of complex engineering systems is decomposed into a set of subproblems in ATC approach, where the desired design specifications or targets are established at the top level and cascaded down to each lower level subproblem. The target transformation is a decomposition process from up to down in a level-by-level way. Moreover, subproblems need to get the responses to match the cascaded targets. If targets cascaded to subsystems are not achievable or compatible, the feedback from subsystems is required, resulting in an iterative target cascading process. Hence, targets and responses are updated and coordinated iteratively to achieve global consistency (see Figure1) [10]. Coordination strategy is needed to implement the communication of target and response values between parents and children problems.Figure 1
Schematic of the ATC process.ATC is a model-based optimization method. There exist two kinds of models, optimal design modelsP and analysis models R. The former is just so-called optimization design problem that includes optimization object, design variables, and constraints. The latter is physical model or mathematical model that is the foundation of the optimization design. Figure 2 shows the data flow between subproblems (parent and children) [19].Figure 2
Information flow between children subproblem and parent subproblem.In Figure2, there exist consistency constraints for each subproblem. In order to make each subproblem independent, consistency constraint is introduced as a kind of equality constraints. For the subproblem at a given level, consistency constraints are as follows:(1)
c
1
′
=
t
i
j
-
r
i
j
=
0
,
c
1
′
=
R
i
j
i
-
1
-
R
i
j
i
=
0
,
c
1
′
=
(
y
i
p
i
-
1
-
y
i
j
i
)
=
0
,
c
2
′
=
t
i
+
1
k
-
r
i
+
1
k
=
0
,
c
2
′
=
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
=
0
,
c
2
′
=
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
=
0
.The general mathematical notation for subproblemsP
i
j in the hierarchy is then(2)
P
i
j
:
min
⁡
f
i
j
w
.
r
.
t
.
k
x
i
j
,
y
i
j
i
,
R
(
i
+
1
)
j
i
,
y
(
i
+
1
)
j
i
where
R
i
j
=
r
i
j
(
R
(
i
+
1
)
j
i
,
x
i
j
,
y
i
j
i
)
s
.
t
.
k
c
1
=
(
R
i
j
i
-
1
-
R
i
j
i
)
=
0
,
c
2
=
(
y
i
p
i
-
1
-
y
i
j
i
)
=
0
,
k
k
c
3
=
(
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
)
=
0
,
k
k
c
4
=
(
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
)
=
0
,
k
k
g
i
j
(
R
i
+
1
k
1
i
,
…
,
R
(
i
+
1
)
k
c
i
j
i
,
x
i
j
,
y
i
j
i
)
≤
0
k
k
h
i
j
(
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
)
=
0
,where P
i
j is the problem formulation of element j at level i. p is used to designate the parent of element j. k designates children of element j. x
i
j is local variables for subproblem P
i
j. r
i
j is the function that calculates responses for P
i
j. R
i
j
i is response of P
i
j. R
i
j
i
-
1 is the target for P
i
j cascaded from the (i
-
1) level parent element. R
(
i
+
1
)
k
i
+
1 is response of the (i
+
1) level children element. R
(
i
+
1
)
j
i
=
(
R
(
i
+
1
)
k
1
i
,
…
,
R
(
i
+
1
)
k
c
i
j
i
) is the target for the (i
+
1) level children element. y
i
j
i is response of linking variable of P
i
j. y
i
p
i
-
1 is the target of linking variable for P
i
j from the (i
-
1) level parent element. y
(
i
+
1
)
k
i
+
1 is response of linking variable of the (i
+
1) level children element. y
(
i
+
1
)
j
i is the target of linking variable for the (i
+
1) level children element from P
i
j. f
i
j is the local objective for P
i
j. g
i
j is function of inequality constraints for P
i
j in negative null form. h
i
j is function of equality constraints for P
i
j in null form.The coordination strategy in ATC deals with the definition of deviation between targets and responses (namely, consistency constraint) and updating weight coefficients therein. By setting appropriate norm and weight coefficients to coordinate each subproblem, the optimal design point can be obtained. The coordination strategy is the key research in ATC, and an appropriate coordination strategy can improve the computation efficiency of ATC. In recent years, some coordination strategies in ATC have been studied.


### 2.1.1. The Coordination Strategy Based on ALF

According to the ALF, the deviation between targets and responses is presented by augmented Lagrangian function, and weight coefficients are presented by the augmented Lagrange multiplier. The coordination way based on ALF is shown as(3)
min
⁡
f
i
j
+
λ
i
j
R
(
R
i
j
i
-
R
i
j
i
-
1
)
+
w
i
j
R
∘
R
i
j
i
-
R
i
j
i
-
1
2
+
λ
i
j
y
(
y
i
p
i
-
1
-
y
i
j
i
)
+
w
i
j
y
∘
y
i
p
i
-
1
-
y
i
j
i
2
+
λ
i
+
1
j
R
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
+
w
(
i
+
1
)
j
R
∘
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
2
+
λ
i
+
1
j
y
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
+
w
i
+
1
j
y
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
2
,
s
.
t
.
g
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
≤
0
k
h
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
=
0
.Weight coefficientsλ
i
j
R, λ
i
j
y, w
i
j
R, and w
i
j
y in (3) are updated in accordance with the following equation:(4)
λ
i
j
R
k
+
1
=
λ
i
j
R
k
+
2
w
i
j
R
k
R
i
j
i
k
-
R
i
j
i
-
1
k
,
λ
i
j
y
k
+
1
=
λ
i
j
y
k
+
2
w
i
j
y
k
y
i
j
i
k
-
y
i
j
i
-
1
k
,
w
i
j
R
k
+
1
=
β
w
i
j
R
k
,
w
i
j
y
k
+
1
=
β
w
i
j
y
k
,where the superscript k indicates iterations and β is step length. Usually, β
≥
1, but 2
<
β
<
3 can accelerate convergence.


### 2.1.2. The Coordination Strategy Based on DQA

DQA is realized by approximating the deviation between targets and responses on the basis of ALF. The specific approximation is(5)
R
i
j
i
-
R
i
j
i
-
1
2
=
R
i
j
i
2
-
2
R
i
j
i
R
i
j
i
-
1
+
R
i
j
i
-
1
2
where
,
R
i
j
i
R
i
j
i
-
1
≅
R
i
j
i
k
-
1
R
i
j
i
-
1
k
-
1
+
R
i
j
i
k
-
1
R
i
j
i
-
R
i
j
i
-
1
k
-
1
+
R
i
j
i
-
1
k
-
1
R
i
j
i
-
R
i
j
i
k
-
1
so
,
R
i
j
i
-
R
i
j
i
-
1
2
=
R
i
j
i
-
1
k
-
1
-
R
i
j
i
2
+
R
i
j
i
-
R
i
j
i
-
1
k
-
1
2
+
R
i
j
i
k
-
1
R
i
j
i
-
1
k
-
1
(
constant
value
)
.Based on DQA, the cross term among the deviation is approximated; that is to say, an order Taylor type expansion is conducted in the last design point. Among them, the weight coefficient is updated according to (5) too.


### 2.1.3. The Coordination Strategy Based on LDF

According to the LDF, the deviation is dealt with by Lagrangian duality function, and weight coefficients are presented by the dual Lagrange multiplier. The coordination of LDF is shown as(6)
min
⁡
f
i
j
+
λ
i
j
R
R
i
j
i
-
R
i
j
i
-
1
+
μ
i
j
R
∘
R
i
j
i
-
R
i
j
i
-
1
2
kkkkkkk
+
λ
i
j
y
y
i
p
i
-
1
-
y
i
j
i
+
μ
i
j
y
∘
y
i
p
i
-
1
-
y
i
j
i
2
s
.
t
.
g
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
≤
0
k
k
k
k
l
h
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
=
0
.Weight coefficientsλ
i
j
R, λ
i
j
y, μ
i
j
R, and μ
i
j
y in (6) are updated in accordance with the following equation:(7)
λ
i
j
R
k
+
1
=
λ
i
j
R
k
+
α
i
j
R
k
R
i
j
i
k
-
R
i
j
i
-
1
k
,
λ
i
j
y
k
+
1
=
λ
i
j
y
k
+
α
i
j
y
k
y
i
j
i
k
-
y
i
j
i
-
1
k
,
μ
i
j
R
k
+
1
=
λ
i
j
R
k
,
μ
i
j
y
k
+
1
=
λ
i
j
y
k
,
α
i
j
R
k
+
1
=
1
+
m
k
+
m
·
1
R
i
j
i
k
-
R
i
j
i
-
1
k
,
α
i
j
y
k
+
1
=
1
+
m
k
+
m
·
1
y
i
j
i
k
-
y
i
j
i
-
1
k
,where the superscript k indicates iterations and m is a positive integers.The difference of these coordination strategies is the way of dealing with consistency constraints and weight coefficients. There exist multiweight coefficients in the above three coordination strategies to increase computational cost. For this reason, the coordination strategy based on quadratic exterior penalty function among which only one weight coefficient exists is put forward, which simplifies the expression of optimization problem greatly, thus improving the efficiency to make ATC have wider application. The coordination strategy based on the quadratic exterior penalty function will be described briefly.


## 2.1.1. The Coordination Strategy Based on ALF

According to the ALF, the deviation between targets and responses is presented by augmented Lagrangian function, and weight coefficients are presented by the augmented Lagrange multiplier. The coordination way based on ALF is shown as(3)
min
⁡
f
i
j
+
λ
i
j
R
(
R
i
j
i
-
R
i
j
i
-
1
)
+
w
i
j
R
∘
R
i
j
i
-
R
i
j
i
-
1
2
+
λ
i
j
y
(
y
i
p
i
-
1
-
y
i
j
i
)
+
w
i
j
y
∘
y
i
p
i
-
1
-
y
i
j
i
2
+
λ
i
+
1
j
R
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
+
w
(
i
+
1
)
j
R
∘
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
2
+
λ
i
+
1
j
y
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
+
w
i
+
1
j
y
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
2
,
s
.
t
.
g
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
≤
0
k
h
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
=
0
.Weight coefficientsλ
i
j
R, λ
i
j
y, w
i
j
R, and w
i
j
y in (3) are updated in accordance with the following equation:(4)
λ
i
j
R
k
+
1
=
λ
i
j
R
k
+
2
w
i
j
R
k
R
i
j
i
k
-
R
i
j
i
-
1
k
,
λ
i
j
y
k
+
1
=
λ
i
j
y
k
+
2
w
i
j
y
k
y
i
j
i
k
-
y
i
j
i
-
1
k
,
w
i
j
R
k
+
1
=
β
w
i
j
R
k
,
w
i
j
y
k
+
1
=
β
w
i
j
y
k
,where the superscript k indicates iterations and β is step length. Usually, β
≥
1, but 2
<
β
<
3 can accelerate convergence.


## 2.1.2. The Coordination Strategy Based on DQA

DQA is realized by approximating the deviation between targets and responses on the basis of ALF. The specific approximation is(5)
R
i
j
i
-
R
i
j
i
-
1
2
=
R
i
j
i
2
-
2
R
i
j
i
R
i
j
i
-
1
+
R
i
j
i
-
1
2
where
,
R
i
j
i
R
i
j
i
-
1
≅
R
i
j
i
k
-
1
R
i
j
i
-
1
k
-
1
+
R
i
j
i
k
-
1
R
i
j
i
-
R
i
j
i
-
1
k
-
1
+
R
i
j
i
-
1
k
-
1
R
i
j
i
-
R
i
j
i
k
-
1
so
,
R
i
j
i
-
R
i
j
i
-
1
2
=
R
i
j
i
-
1
k
-
1
-
R
i
j
i
2
+
R
i
j
i
-
R
i
j
i
-
1
k
-
1
2
+
R
i
j
i
k
-
1
R
i
j
i
-
1
k
-
1
(
constant
value
)
.Based on DQA, the cross term among the deviation is approximated; that is to say, an order Taylor type expansion is conducted in the last design point. Among them, the weight coefficient is updated according to (5) too.


## 2.1.3. The Coordination Strategy Based on LDF

According to the LDF, the deviation is dealt with by Lagrangian duality function, and weight coefficients are presented by the dual Lagrange multiplier. The coordination of LDF is shown as(6)
min
⁡
f
i
j
+
λ
i
j
R
R
i
j
i
-
R
i
j
i
-
1
+
μ
i
j
R
∘
R
i
j
i
-
R
i
j
i
-
1
2
kkkkkkk
+
λ
i
j
y
y
i
p
i
-
1
-
y
i
j
i
+
μ
i
j
y
∘
y
i
p
i
-
1
-
y
i
j
i
2
s
.
t
.
g
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
≤
0
k
k
k
k
l
h
i
j
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
=
0
.Weight coefficientsλ
i
j
R, λ
i
j
y, μ
i
j
R, and μ
i
j
y in (6) are updated in accordance with the following equation:(7)
λ
i
j
R
k
+
1
=
λ
i
j
R
k
+
α
i
j
R
k
R
i
j
i
k
-
R
i
j
i
-
1
k
,
λ
i
j
y
k
+
1
=
λ
i
j
y
k
+
α
i
j
y
k
y
i
j
i
k
-
y
i
j
i
-
1
k
,
μ
i
j
R
k
+
1
=
λ
i
j
R
k
,
μ
i
j
y
k
+
1
=
λ
i
j
y
k
,
α
i
j
R
k
+
1
=
1
+
m
k
+
m
·
1
R
i
j
i
k
-
R
i
j
i
-
1
k
,
α
i
j
y
k
+
1
=
1
+
m
k
+
m
·
1
y
i
j
i
k
-
y
i
j
i
-
1
k
,where the superscript k indicates iterations and m is a positive integers.The difference of these coordination strategies is the way of dealing with consistency constraints and weight coefficients. There exist multiweight coefficients in the above three coordination strategies to increase computational cost. For this reason, the coordination strategy based on quadratic exterior penalty function among which only one weight coefficient exists is put forward, which simplifies the expression of optimization problem greatly, thus improving the efficiency to make ATC have wider application. The coordination strategy based on the quadratic exterior penalty function will be described briefly.


## 2.2. The QEPF-ATC Approach

A general equality-constrained optimization problem can be defined as follows:(8)
min
⁡
f
(
x
)
s
.
t
.
lk
c
i
(
x
)
=
0
,
i
=
1,2
,
…
,
m
k
k
k
k
l
k
x
lb
≤
x
≤
x
up
,where x is the vector of the design variables with x
lb and x
up as its lower and upper bounds, respectively. The function f
(
x
) represents the objective function. c
i
(
x
) is the equality constraints.According to QEPF, consistency constraints are added into the objective function by introducing a penalty factor, and the penalty factor is used as the weighting coefficient to realize the coordination. The quadratic exterior penalty function is defined as(9)
F
x
,
ρ
=
f
x
+
ρ
2
∑
i
=
1
m
c
i
x
2
=
f
x
+
ρ
2
c
x
T
c
x
,where the weight coefficient (penalty factor) is updated in accordance with (10). The value of λ influences the convergence efficiency and computational efficiency so that it is necessary to set up appropriate λ to obtain much better solution efficiency. Consider (10)
ρ
k
+
1
=
λ
ρ
k
,
λ
=
1.4

~


10
.The mathematical formulation of QEPF-ATC is(11)
P
i
j
:
min
⁡
f
i
j
+
ρ
2
R
i
j
i
-
1
-
R
i
j
i
2
+
y
i
p
i
-
1
-
y
i
j
i
2
k
+
R
i
+
1
j
i
-
R
i
+
1
k
i
+
1
2
k
+
y
i
+
1
j
i
-
y
i
+
1
k
i
+
1
2
s
.
t
.
g
i
j
≤
0
,
h
i
j
=
0
l
R
i
j
=
r
i
j
x
i
j
¯
,
l
x
i
j
¯
=
R
i
+
1
k
1
i
,
…
,
R
i
+
1
k
c
i
j
i
,
x
i
j
,
y
i
j
i
.Generally, the much bigger the penalty factor is the much closer the optimization results are to the optimal solution of the original problem. But, in the practical application, the value of the penalty factor must be moderate when considering the computational cost, due to the fact that if the penalty factor is too small, the optimization results may get bigger error from the original problem to bring low computational efficiency. Otherwise, too large penalty factor would increase the computational expense.Figure3 shows the solving process of the QEPF-ATC approach, and the main steps are as follows:(1)
decompose the original problem and confirm design variables and constraints of system and subsystems;
(2)
initialize design variables (including local design variables and linking variables) andρ
(
0
) and transfer the corresponding initial value to the subsystem as the target;
(3)
solve each subsystem and pass back current results to system level;
(4)
optimize the system problem on the basis of the third step and transfer new target to each subsystem;
(5)
obtain the optimization result if the convergence condition is met; otherwise, update the weight coefficientρ according to (10) and go back to the third step until the convergence condition is met.Figure 3
The solving process of QEPF-ATC approach.


## 3. The Proposed QEPF-ATC and Kriging Model Combined Approach

### 3.1. Kriging Model

Kriging models [15, 16] have their origins in mining and geostatistical applications involving spatially and temporally correlated data. Kriging models combine a global model with localized departures:(12)
y
x
=
f
x
+
z
x
,where y
(
x
) is the unknown function of interest, f
(
x
) is the known approximation (usually polynomial) function, and z
(
x
) is the realization of a stochastic process with mean zero, variance σ
2, and nonzero covariance. The f
(
x
) term is similar to a polynomial response surface, providing a “global” model of the design space. The covariance of z
(
x
) is given by(13)
cov
z
x
i
,
z
x
j
=
σ
2
R
x
i
,
x
j
,where R is the correlation. R
(
x
i
,
x
j
) is the correlation function between any two of the n sampled data points x
i and x
j, and it can be calculated by(14)
R
x
i
,
x
j
=
exp
⁡
⁡
-
∑
k
=
1
n
θ
k
x
i
k
-
x
j
k
2
.Generally,f
(
x
) is treated as a constant β, then predicted estimates, y
⌢
(
x
), of the response y
(
x
) at untried values of x are given by(15)
y
⌢
x
=
β
⌢
+
r
T
x
R
-
1
Y
-
F
β
⌢
,where y is the column vector of length n that contains the sample values of the response, and f is a column vector of length n that is filled with ones when f
(
x
) is taken as a constant. In (10), r
T
(
x
) is the correlation vector of length n between an untried x and the sampled data points.Consider(16)
r
T
x
=
R
x
,
x
1
,
R
x
,
x
2
,
R
x
,
x
3
,
…
,
R
x
,
x
n
T
.In (10), β
⌢ is estimated using (12):(17)
β
⌢
=
F
T
R
-
1
F
-
1
F
T
R
-
1
YInR, θ
k
=
max
⁡
⁡
{
-
[
n
ln
⁡
(
σ
2
)
+
ln
⁡
R
]
/
2
},(18)
σ
2
=
Y
-
F
β
⌢
T
R
-
1
(
Y
-
F
β
⌢
)
n
.


### 3.2. The QEPF-ATC and Kriging Model Combined Approach

In QEPF-ATC, the deviation between the target and response is processed in the form of sum of squares, as shown in (11). The sum of squares function is usually not smooth and nonlinear. And discrepancy of variables in the neighborhood of the current solution will lead to fluctuation of the deviation, which will increase the number of iterations and so forth. Kriging model is appropriate for high nonlinear problem with better comprehensive performance, high approximation accuracy, and high robustness [20–22]. Moreover, fewer sample points also can get higher precision approximation model. So Kriging model is adopted in this paper.The solving process of the QEPF-ATC and Kriging model combined approach is shown in Figure4, and the main steps are as follows:(1)
decompose the original problem and confirm design variables and constraints of system and subsystems;
(2)
select a set of sample points using the Latin hypercube sampling (LHS) method;
(3)
analyse each subsystem to obtain responses at all sampling points of the subsystems;
(4)
build Kriging model based on the sample points and the corresponding responses of each subsystem;
(5)
evaluate the Kriging model by error analysis; if the precision conditions are not met, go back to step(
2
) to build a new Kriging model by adding some new sample points;
(6)
initialize variables and transfer the corresponding initial value to the subsystem as the target and optimize each subproblem based on the created Kriging model;
(7)
pass back current results to system level, and implement the system level optimization;
(8)
obtain results if the convergence condition is met; otherwise, go to step(
6
) until the convergence condition is met.Figure 4
The solving process of the QEPF-ATC and Kriging model combined approach.


## 3.1. Kriging Model

Kriging models [15, 16] have their origins in mining and geostatistical applications involving spatially and temporally correlated data. Kriging models combine a global model with localized departures:(12)
y
x
=
f
x
+
z
x
,where y
(
x
) is the unknown function of interest, f
(
x
) is the known approximation (usually polynomial) function, and z
(
x
) is the realization of a stochastic process with mean zero, variance σ
2, and nonzero covariance. The f
(
x
) term is similar to a polynomial response surface, providing a “global” model of the design space. The covariance of z
(
x
) is given by(13)
cov
z
x
i
,
z
x
j
=
σ
2
R
x
i
,
x
j
,where R is the correlation. R
(
x
i
,
x
j
) is the correlation function between any two of the n sampled data points x
i and x
j, and it can be calculated by(14)
R
x
i
,
x
j
=
exp
⁡
⁡
-
∑
k
=
1
n
θ
k
x
i
k
-
x
j
k
2
.Generally,f
(
x
) is treated as a constant β, then predicted estimates, y
⌢
(
x
), of the response y
(
x
) at untried values of x are given by(15)
y
⌢
x
=
β
⌢
+
r
T
x
R
-
1
Y
-
F
β
⌢
,where y is the column vector of length n that contains the sample values of the response, and f is a column vector of length n that is filled with ones when f
(
x
) is taken as a constant. In (10), r
T
(
x
) is the correlation vector of length n between an untried x and the sampled data points.Consider(16)
r
T
x
=
R
x
,
x
1
,
R
x
,
x
2
,
R
x
,
x
3
,
…
,
R
x
,
x
n
T
.In (10), β
⌢ is estimated using (12):(17)
β
⌢
=
F
T
R
-
1
F
-
1
F
T
R
-
1
YInR, θ
k
=
max
⁡
⁡
{
-
[
n
ln
⁡
(
σ
2
)
+
ln
⁡
R
]
/
2
},(18)
σ
2
=
Y
-
F
β
⌢
T
R
-
1
(
Y
-
F
β
⌢
)
n
.


## 3.2. The QEPF-ATC and Kriging Model Combined Approach

In QEPF-ATC, the deviation between the target and response is processed in the form of sum of squares, as shown in (11). The sum of squares function is usually not smooth and nonlinear. And discrepancy of variables in the neighborhood of the current solution will lead to fluctuation of the deviation, which will increase the number of iterations and so forth. Kriging model is appropriate for high nonlinear problem with better comprehensive performance, high approximation accuracy, and high robustness [20–22]. Moreover, fewer sample points also can get higher precision approximation model. So Kriging model is adopted in this paper.The solving process of the QEPF-ATC and Kriging model combined approach is shown in Figure4, and the main steps are as follows:(1)
decompose the original problem and confirm design variables and constraints of system and subsystems;
(2)
select a set of sample points using the Latin hypercube sampling (LHS) method;
(3)
analyse each subsystem to obtain responses at all sampling points of the subsystems;
(4)
build Kriging model based on the sample points and the corresponding responses of each subsystem;
(5)
evaluate the Kriging model by error analysis; if the precision conditions are not met, go back to step(
2
) to build a new Kriging model by adding some new sample points;
(6)
initialize variables and transfer the corresponding initial value to the subsystem as the target and optimize each subproblem based on the created Kriging model;
(7)
pass back current results to system level, and implement the system level optimization;
(8)
obtain results if the convergence condition is met; otherwise, go to step(
6
) until the convergence condition is met.Figure 4
The solving process of the QEPF-ATC and Kriging model combined approach.


## 4. Case Studies

### 4.1. The Geometric Programming

The geometric programming [23] that has only the global optimal solution is employed to verify the usability of QEPF-ATC. Its mathematical expression is shown as(19)
min
⁡
f
=
x
1
2
+
x
2
2
,
s
.
t
.
k
g
1
:
x
3
-
2
+
x
4
2
x
5
2
≤
1
,
g
2
:
x
5
2
+
x
6
-
2
x
7
2
≤
1
,
k
g
3
:
x
8
2
+
x
9
2
x
11
2
≤
1
,
k
g
4
:
x
8
-
2
+
x
10
2
x
11
2
≤
1
,
g
5
:
x
11
2
+
x
12
-
2
x
13
2
≤
1
,
k
g
6
:
x
11
2
+
x
12
2
x
14
2
≤
1
,
k
h
1
=
x
1
2
-
x
3
2
-
x
4
-
2
-
x
5
2
=
0
,
k
h
2
=
x
2
2
-
x
5
2
-
x
6
2
-
x
7
2
=
0
,
k
h
3
=
x
3
2
-
x
8
2
-
x
9
-
2
-
x
10
-
2
-
x
11
2
=
0
,
k
h
4
=
x
6
2
-
x
11
2
-
x
12
2
-
x
13
2
-
x
14
2
=
0
.Figure5 shows the solution framework of the QEPF-ATC approach for the problem, where x
11 is the linking variable (coupling variable) coordinated by the system level. x
11
B
U, x
11
C
U, x
3
U, and x
6
U are the targets of each subproblem passed from the system level. x
11
B
L and x
11
C
L are the responses of the linking variable coming from the subproblems B and C. x
3
L and x
6
L are the responses of B and C.Figure 5
The solving framework of QEPF-ATC for geometric programming.Solve the geometric programming problem according to QEPF-ATC. Initialize the design pointx
0
=
1,1
,
1,1
,
1,1
,
1,1
,
1,1
,
1,1
,
1,1
T, and take λ
=
10. Here, the acceptable inconsistency tolerance is 0.01 for every response variable and linking variable. Table 1 shows the optimal results based on QEPF-ATC, three coordination strategies (ALF, DQA, and LDF), and AAO (All-At-Once).Table 1
The optimal results.
QEPF-ATC
ALF
DQA
LDF
AAO
x
1
2.9049
2.85
2.86
2.87
2.84
x
2
3.0551
3.04
3.10
3.05
3.09
x
3
2.4930
2.48
2.49
2.43
2.36
x
4
0.8805
0.86
0.82
0.79
0.76
x
5
0.9676
0.94
0.95
0.89
0.87
x
6
2.7103
2.72
2.75
2.78
2.81
x
7
1.0357
1.02
0.95
1.03
0.94
x
8
0.9830
0.98
0.98
1.01
0.97
x
9
0.7572
0.78
0.79
0.82
0.87
x
10
0.7104
0.76
0.79
0.81
0.80
x
11
1.2408
1.27
1.26
1.25
1.30
x
12
0.8409
0.86
0.85
0.83
0.84
x
13
1.7088
1.73
1.74
1.71
1.76
x
14
1.4875
1.52
1.49
1.50
1.55
f
17.7721
17.36
17.78
17.53
17.59
Runtime (s)
289
355
336
368According to Table1, the results of the QEPF-ATC approach are very close to the optimal value of the original problem, and the runtime is less than the other three methods.λ affects the convergence efficiency and computational efficiency of the QEPF-ATC approach. In order to obtain more suitable λ to improve the efficiency of the QEPF-ATC approach, different values of λ are selected to solve the problem. Table 2 shows the results when different values of λ are selected.Table 2
The optimal results based onλ.
λ
(x
1, x
2, x
3, x
4, x
5, x
6, x
7, x
11)
(x
3, x
8, x
9, x
10, x
11)
(x
6, x
11, x
12, x
13, x
14)
f
Iterations
1.4
(2.7556, 3.0804, 2.2975, 0.8131, 0.9171, 2.7718, 0.9840, 1.2678)
(2.3012, 0.9688, 0.8427, 0.7340,  1.2682)
(2.7638, 1.2732, 0.7447, 1.8411, 1.4655)
17.5992
78
4
(2.9716, 2.9567, 2.5227, 0.7598, 0.8570, 2.6708, 0.9352, 1.2343)
(2.5230, 0.9911, 0.7320, 0.7072, 1.2321)
(2.6791, 1.2298, 0.8409, 1.6923, 1.4686)
17.5728
27
7
(2.9662, 2.9757, 2.5012, 0.7204, 0.8231, 2.7169, 0.9014, 1.2386)
(2.5011, 0.9863, 0.7492, 0.7115, 1.2386)
(2.7168, 1.2385, 0.8408, 1.7122, 1.4915)
17.6531
17
10
(2.9049, 3.0551, 2.4930, 0.8805, 0.9676, 2.7103, 1.0357, 1.2408)
(2.4930, 0.9830, 0.7572, 0.7104, 1.2408)
(2.7102, 1.2409, 0.8409, 1.7088, 1.4875)
17.7721
5In Table2, the QEPF-ATC approach implement fast convergence when λ is set to different values, and the optimal results are close to the results of the initial problems. With the increase of λ, the number of iterations decreases and the solution error gets bigger (excessive weight). Generally, 1.4
≤
λ
≤
10 is strictly necessary, but typically, 4
<
λ
<
7 is recommended to speed up convergence and reduce computational cost and get better optimal results.


### 4.2. The Speed Reducer Design

The speed reducer design case represents the design of a simple gearbox and is a multidisciplinary problem comprising the coupling between gear design and shaft design. This has been used as a testing problem for MDO method in the literature from NASA Langley research center [24]. The multidisciplinary systems and notations of the reducer are shown in Figure 6, and the modelling is shown in Figure 7. The design objective is to minimize the speed reducer weight while satisfying a number of constraints posed by gear and shaft parts. There are seven variables in the speed reducer design optimization. Table 3 lists a brief description of the variables and gives their reference values.Table 3
Design variables for speed reducer.
Variables
Range
Units
x
1
Gear face width
2.6
≤
x
1
≤
3.6
cm
x
2
Teeth module
0.7
≤
x
2
≤
0.8
cm
x
3
Number of teeth of pinion (integer variable)
17
≤
x
3
≤
28
/
x
4
Distance between bearings 1
7.3
≤
x
4
≤
8.3
cm
x
5
Distance between bearings 2
7.3
≤
x
5
≤
8.3
cm
x
6
Diameter of shaft 1
2.9
≤
x
6
≤
3.9
cm
x
7
Diameter of shaft 1
5.0
≤
x
7
≤
5.5
cm
x
1, x
2, and x
3 are the design variables of the gear subproblem. x
4 and x
5 are the design variables of the shaft subproblem. x
6 and x
7 are linking variables (coupling variables).Figure 6
The multidisciplinary systems and notations of reducer.Figure 7
The decomposition of speed reducer problem.According to the decomposition method of ATC, where a system is partitioned by object, the speed reducer problem is decomposed into gear subproblem and shaft subproblem (shown in Figure7). The design optimization problem is defined as follows:(20)
min
⁡
f
=
f
1
x
+
f
2
x
,
k
l
f
1
(
x
)
=
0.7854
x
1
x
2
2
·
(
3.3333
x
3
2
+
14.9334
x
3
-
43.0934
)
-
1.5079
x
1
x
6
2
+
x
7
2
,
k
l
f
2
(
x
)
=
7.477
(
x
6
3
+
x
7
3
)
+
0.7854
(
x
4
x
6
2
+
x
5
x
7
2
)
s
.
t
.
k
g
1
=
27
x
1
x
2
2
x
3
-
1
≤
0
,
g
2
=
397.5
x
1
x
2
2
x
3
2
-
1
≤
0
,
k
l
g
3
=
1.93
x
4
3
x
2
x
3
x
6
4
-
1
≤
0
,
g
4
=
1.93
x
5
3
x
2
x
3
x
7
4
-
1
≤
0
,
k
l
g
5
=
10
x
6
-
3
745
x
2
-
1
x
3
-
1
x
4
+
1.69
×
1
0
7
kkkkkkkkkk
-
1100
≤
0
,
k
l
g
6
=
10
x
7
-
3
745
x
2
-
1
x
3
-
1
x
5
+
1.575
×
1
0
8
kkkkkkkkk
-
850
≤
0
,
k
l
g
7
=
1.5
x
6
+
1.9
x
4
-
1
≤
0
,
k
l
g
8
=
1.5
x
7
+
1.9
x
5
-
1
≤
0
,
k
l
g
9
=
x
2
x
3
-
40
≤
0
,
g
10
=
5
-
x
1
x
2
≤
0
,
k
l
g
11
=
x
1
x
2
-
12
≤
0
,where f
1 is the gear weight and f
2 is the shaft weight. And g
1 and g
2, are respectively, gear bending stress and contact stress constraints, g
3 and g
4 are, respectively, torsional distortion constraints of shaft 1 and shaft 2, g
5 and g
6 are stress constraints of the shaft 1 and shaft 2, g
7 and g
8 are experience constraints of shaft 1 and shaft 2, and g
9, g
10, and g
11 are gear geometry constraints.x
1, x
2, and x
3 are the design variables of the gear subproblem. x
4 and x
5 are the design variables of the shaft subproblem. x
6 and x
7 are linking variables (coupling variables). According to the QEPF-ATC and Kriging model combined approach, three variables x
8, x
9, and x
10 are introduced and meet the following formula:(21)
h
0
=
x
8
+
x
9
-
x
10
=
0
,
h
1
=
x
8
-
f
1
x
=
0
,
h
2
=
x
9
-
f
2
x
=
0
,where h
0 represents the analysis model of system level and h
1 and h
2 represent the analysis models of subproblem 1 and subproblem 2, respectively.The solving framework of the QEPF-ATC and Kriging model combined approach for gear reducer is shown in Figure8.Figure 8
The solving framework of the QEPF-ATC and Kriging model combined approach for the gear reducer.Creating the Kriging models of the shaft subsystem is introduced in detail. Firstly, select 45 sample points by LHS method to create the Kriging models to conduct the subsystem analysis, and the created models are shown in Figure9.Figure 9
Kriging models of the shaft subsystem.Then, select 20 sample points by random sampling method to evaluate the created Kriging models. Error analysis for Kriging models of the shaft subsystem is shown in Figure10, and it can be seen that the predicted values obtained by Kriging models are basically equal to the actual values.Figure 10
Error analysis for Kriging models of the shaft subsystem.There are four credibility evaluation standards for approximate models, that is, the average, the maximum, the root mean square, and theR-squared. For the average, the maximum, and the root mean square, the much smaller the values are the much better the credibility of approximate models is. The default threshold of the average and the root mean square is 0.2, and the maximum is 0.3. For the R-squared, the much bigger the value is, the much better the credibility of approximate models is, and the default threshold is 0.9. Figure 11 describes the values of the four credibility evaluation standards for the created Kriging models of the shaft subsystem, and it can be seen that the four credibility evaluation standards are all met. So the created Kriging models are reliable enough and can be used to replace the original models with high accuracy in the next solving process.Figure 11
The error data of Kriging models of subsystem 2.Five initial points are selected and the sequential quadratic programming (SQP) method is adopted to solve the gear reducer problem, and the results provided by the QEPF-ATC and Kriging model combined approach and the QEPF-ATC approach are shown in Tables4 and 5, respectively. It can be seen that results of the QEPF-ATC and Kriging model combined approach and the QEPF-ATC approach are very close, and the best result provided by the QEPF-ATC and Kriging model combined approach is closer to the optimal value of the gear reducer problem compared to that provided by the QEPF-ATC approach. In terms of the iterations and runtime, the QEPF-ATC and Kriging model combined approach saves computing time by more than one half compared with QEPF-ATC approach.Table 4
The optimal results of the QEPF-ATC and Kriging model combined approach.
Initial point
f
min
⁡
Iterations
Runtime (s)
(x
1, x
2, x
3, x
4, x
5, x
6, x
7)
1
2994.53432
26
138
(3.504, 0.698, 17, 7.300, 7.845, 3.410, 5.283)
2
2994.53217
25
124
(3.501, 0.700, 17, 7.310, 7.763, 3.353, 5.286)
3
2995.21246
18
113
(3.498, 0.699, 17, 7.306, 7.755, 3.312, 5.279)
4
2995.64235
10
36
(3.499, 0.702, 17, 7.297, 7.718, 3.350, 5.282)
5
2994.33263
24
121
(3.502, 0.700, 17, 7.298, 7.720, 3.347, 5.287)Table 5
The optimal results of the QEPF-ATC approach.
Initial point
f
min
⁡
Iterations
Runtime (s)
(x
1, x
2, x
3, x
4, x
5, x
6, x
7)
1
2994.63432
42
272
(3.500, 0.699, 17, 7.300, 7.825, 3.392, 5.282)
2
2995.51327
38
215
(3.500, 0.703, 17, 7.305, 7.746, 3.353, 5.287)
3
2996.12356
30
198
(3.497, 0.696, 17, 7.300, 7.713, 3.326, 5.275)
4
2995.75363
17
97
(3.499, 0.700, 17, 7.299, 7.715, 3.349, 5.282)
5
2995.21554
52
306
(3.500, 0.700, 17, 7.300, 7.719, 3.350, 5.286)


## 4.1. The Geometric Programming

The geometric programming [23] that has only the global optimal solution is employed to verify the usability of QEPF-ATC. Its mathematical expression is shown as(19)
min
⁡
f
=
x
1
2
+
x
2
2
,
s
.
t
.
k
g
1
:
x
3
-
2
+
x
4
2
x
5
2
≤
1
,
g
2
:
x
5
2
+
x
6
-
2
x
7
2
≤
1
,
k
g
3
:
x
8
2
+
x
9
2
x
11
2
≤
1
,
k
g
4
:
x
8
-
2
+
x
10
2
x
11
2
≤
1
,
g
5
:
x
11
2
+
x
12
-
2
x
13
2
≤
1
,
k
g
6
:
x
11
2
+
x
12
2
x
14
2
≤
1
,
k
h
1
=
x
1
2
-
x
3
2
-
x
4
-
2
-
x
5
2
=
0
,
k
h
2
=
x
2
2
-
x
5
2
-
x
6
2
-
x
7
2
=
0
,
k
h
3
=
x
3
2
-
x
8
2
-
x
9
-
2
-
x
10
-
2
-
x
11
2
=
0
,
k
h
4
=
x
6
2
-
x
11
2
-
x
12
2
-
x
13
2
-
x
14
2
=
0
.Figure5 shows the solution framework of the QEPF-ATC approach for the problem, where x
11 is the linking variable (coupling variable) coordinated by the system level. x
11
B
U, x
11
C
U, x
3
U, and x
6
U are the targets of each subproblem passed from the system level. x
11
B
L and x
11
C
L are the responses of the linking variable coming from the subproblems B and C. x
3
L and x
6
L are the responses of B and C.Figure 5
The solving framework of QEPF-ATC for geometric programming.Solve the geometric programming problem according to QEPF-ATC. Initialize the design pointx
0
=
1,1
,
1,1
,
1,1
,
1,1
,
1,1
,
1,1
,
1,1
T, and take λ
=
10. Here, the acceptable inconsistency tolerance is 0.01 for every response variable and linking variable. Table 1 shows the optimal results based on QEPF-ATC, three coordination strategies (ALF, DQA, and LDF), and AAO (All-At-Once).Table 1
The optimal results.
QEPF-ATC
ALF
DQA
LDF
AAO
x
1
2.9049
2.85
2.86
2.87
2.84
x
2
3.0551
3.04
3.10
3.05
3.09
x
3
2.4930
2.48
2.49
2.43
2.36
x
4
0.8805
0.86
0.82
0.79
0.76
x
5
0.9676
0.94
0.95
0.89
0.87
x
6
2.7103
2.72
2.75
2.78
2.81
x
7
1.0357
1.02
0.95
1.03
0.94
x
8
0.9830
0.98
0.98
1.01
0.97
x
9
0.7572
0.78
0.79
0.82
0.87
x
10
0.7104
0.76
0.79
0.81
0.80
x
11
1.2408
1.27
1.26
1.25
1.30
x
12
0.8409
0.86
0.85
0.83
0.84
x
13
1.7088
1.73
1.74
1.71
1.76
x
14
1.4875
1.52
1.49
1.50
1.55
f
17.7721
17.36
17.78
17.53
17.59
Runtime (s)
289
355
336
368According to Table1, the results of the QEPF-ATC approach are very close to the optimal value of the original problem, and the runtime is less than the other three methods.λ affects the convergence efficiency and computational efficiency of the QEPF-ATC approach. In order to obtain more suitable λ to improve the efficiency of the QEPF-ATC approach, different values of λ are selected to solve the problem. Table 2 shows the results when different values of λ are selected.Table 2
The optimal results based onλ.
λ
(x
1, x
2, x
3, x
4, x
5, x
6, x
7, x
11)
(x
3, x
8, x
9, x
10, x
11)
(x
6, x
11, x
12, x
13, x
14)
f
Iterations
1.4
(2.7556, 3.0804, 2.2975, 0.8131, 0.9171, 2.7718, 0.9840, 1.2678)
(2.3012, 0.9688, 0.8427, 0.7340,  1.2682)
(2.7638, 1.2732, 0.7447, 1.8411, 1.4655)
17.5992
78
4
(2.9716, 2.9567, 2.5227, 0.7598, 0.8570, 2.6708, 0.9352, 1.2343)
(2.5230, 0.9911, 0.7320, 0.7072, 1.2321)
(2.6791, 1.2298, 0.8409, 1.6923, 1.4686)
17.5728
27
7
(2.9662, 2.9757, 2.5012, 0.7204, 0.8231, 2.7169, 0.9014, 1.2386)
(2.5011, 0.9863, 0.7492, 0.7115, 1.2386)
(2.7168, 1.2385, 0.8408, 1.7122, 1.4915)
17.6531
17
10
(2.9049, 3.0551, 2.4930, 0.8805, 0.9676, 2.7103, 1.0357, 1.2408)
(2.4930, 0.9830, 0.7572, 0.7104, 1.2408)
(2.7102, 1.2409, 0.8409, 1.7088, 1.4875)
17.7721
5In Table2, the QEPF-ATC approach implement fast convergence when λ is set to different values, and the optimal results are close to the results of the initial problems. With the increase of λ, the number of iterations decreases and the solution error gets bigger (excessive weight). Generally, 1.4
≤
λ
≤
10 is strictly necessary, but typically, 4
<
λ
<
7 is recommended to speed up convergence and reduce computational cost and get better optimal results.


## 4.2. The Speed Reducer Design

The speed reducer design case represents the design of a simple gearbox and is a multidisciplinary problem comprising the coupling between gear design and shaft design. This has been used as a testing problem for MDO method in the literature from NASA Langley research center [24]. The multidisciplinary systems and notations of the reducer are shown in Figure 6, and the modelling is shown in Figure 7. The design objective is to minimize the speed reducer weight while satisfying a number of constraints posed by gear and shaft parts. There are seven variables in the speed reducer design optimization. Table 3 lists a brief description of the variables and gives their reference values.Table 3
Design variables for speed reducer.
Variables
Range
Units
x
1
Gear face width
2.6
≤
x
1
≤
3.6
cm
x
2
Teeth module
0.7
≤
x
2
≤
0.8
cm
x
3
Number of teeth of pinion (integer variable)
17
≤
x
3
≤
28
/
x
4
Distance between bearings 1
7.3
≤
x
4
≤
8.3
cm
x
5
Distance between bearings 2
7.3
≤
x
5
≤
8.3
cm
x
6
Diameter of shaft 1
2.9
≤
x
6
≤
3.9
cm
x
7
Diameter of shaft 1
5.0
≤
x
7
≤
5.5
cm
x
1, x
2, and x
3 are the design variables of the gear subproblem. x
4 and x
5 are the design variables of the shaft subproblem. x
6 and x
7 are linking variables (coupling variables).Figure 6
The multidisciplinary systems and notations of reducer.Figure 7
The decomposition of speed reducer problem.According to the decomposition method of ATC, where a system is partitioned by object, the speed reducer problem is decomposed into gear subproblem and shaft subproblem (shown in Figure7). The design optimization problem is defined as follows:(20)
min
⁡
f
=
f
1
x
+
f
2
x
,
k
l
f
1
(
x
)
=
0.7854
x
1
x
2
2
·
(
3.3333
x
3
2
+
14.9334
x
3
-
43.0934
)
-
1.5079
x
1
x
6
2
+
x
7
2
,
k
l
f
2
(
x
)
=
7.477
(
x
6
3
+
x
7
3
)
+
0.7854
(
x
4
x
6
2
+
x
5
x
7
2
)
s
.
t
.
k
g
1
=
27
x
1
x
2
2
x
3
-
1
≤
0
,
g
2
=
397.5
x
1
x
2
2
x
3
2
-
1
≤
0
,
k
l
g
3
=
1.93
x
4
3
x
2
x
3
x
6
4
-
1
≤
0
,
g
4
=
1.93
x
5
3
x
2
x
3
x
7
4
-
1
≤
0
,
k
l
g
5
=
10
x
6
-
3
745
x
2
-
1
x
3
-
1
x
4
+
1.69
×
1
0
7
kkkkkkkkkk
-
1100
≤
0
,
k
l
g
6
=
10
x
7
-
3
745
x
2
-
1
x
3
-
1
x
5
+
1.575
×
1
0
8
kkkkkkkkk
-
850
≤
0
,
k
l
g
7
=
1.5
x
6
+
1.9
x
4
-
1
≤
0
,
k
l
g
8
=
1.5
x
7
+
1.9
x
5
-
1
≤
0
,
k
l
g
9
=
x
2
x
3
-
40
≤
0
,
g
10
=
5
-
x
1
x
2
≤
0
,
k
l
g
11
=
x
1
x
2
-
12
≤
0
,where f
1 is the gear weight and f
2 is the shaft weight. And g
1 and g
2, are respectively, gear bending stress and contact stress constraints, g
3 and g
4 are, respectively, torsional distortion constraints of shaft 1 and shaft 2, g
5 and g
6 are stress constraints of the shaft 1 and shaft 2, g
7 and g
8 are experience constraints of shaft 1 and shaft 2, and g
9, g
10, and g
11 are gear geometry constraints.x
1, x
2, and x
3 are the design variables of the gear subproblem. x
4 and x
5 are the design variables of the shaft subproblem. x
6 and x
7 are linking variables (coupling variables). According to the QEPF-ATC and Kriging model combined approach, three variables x
8, x
9, and x
10 are introduced and meet the following formula:(21)
h
0
=
x
8
+
x
9
-
x
10
=
0
,
h
1
=
x
8
-
f
1
x
=
0
,
h
2
=
x
9
-
f
2
x
=
0
,where h
0 represents the analysis model of system level and h
1 and h
2 represent the analysis models of subproblem 1 and subproblem 2, respectively.The solving framework of the QEPF-ATC and Kriging model combined approach for gear reducer is shown in Figure8.Figure 8
The solving framework of the QEPF-ATC and Kriging model combined approach for the gear reducer.Creating the Kriging models of the shaft subsystem is introduced in detail. Firstly, select 45 sample points by LHS method to create the Kriging models to conduct the subsystem analysis, and the created models are shown in Figure9.Figure 9
Kriging models of the shaft subsystem.Then, select 20 sample points by random sampling method to evaluate the created Kriging models. Error analysis for Kriging models of the shaft subsystem is shown in Figure10, and it can be seen that the predicted values obtained by Kriging models are basically equal to the actual values.Figure 10
Error analysis for Kriging models of the shaft subsystem.There are four credibility evaluation standards for approximate models, that is, the average, the maximum, the root mean square, and theR-squared. For the average, the maximum, and the root mean square, the much smaller the values are the much better the credibility of approximate models is. The default threshold of the average and the root mean square is 0.2, and the maximum is 0.3. For the R-squared, the much bigger the value is, the much better the credibility of approximate models is, and the default threshold is 0.9. Figure 11 describes the values of the four credibility evaluation standards for the created Kriging models of the shaft subsystem, and it can be seen that the four credibility evaluation standards are all met. So the created Kriging models are reliable enough and can be used to replace the original models with high accuracy in the next solving process.Figure 11
The error data of Kriging models of subsystem 2.Five initial points are selected and the sequential quadratic programming (SQP) method is adopted to solve the gear reducer problem, and the results provided by the QEPF-ATC and Kriging model combined approach and the QEPF-ATC approach are shown in Tables4 and 5, respectively. It can be seen that results of the QEPF-ATC and Kriging model combined approach and the QEPF-ATC approach are very close, and the best result provided by the QEPF-ATC and Kriging model combined approach is closer to the optimal value of the gear reducer problem compared to that provided by the QEPF-ATC approach. In terms of the iterations and runtime, the QEPF-ATC and Kriging model combined approach saves computing time by more than one half compared with QEPF-ATC approach.Table 4
The optimal results of the QEPF-ATC and Kriging model combined approach.
Initial point
f
min
⁡
Iterations
Runtime (s)
(x
1, x
2, x
3, x
4, x
5, x
6, x
7)
1
2994.53432
26
138
(3.504, 0.698, 17, 7.300, 7.845, 3.410, 5.283)
2
2994.53217
25
124
(3.501, 0.700, 17, 7.310, 7.763, 3.353, 5.286)
3
2995.21246
18
113
(3.498, 0.699, 17, 7.306, 7.755, 3.312, 5.279)
4
2995.64235
10
36
(3.499, 0.702, 17, 7.297, 7.718, 3.350, 5.282)
5
2994.33263
24
121
(3.502, 0.700, 17, 7.298, 7.720, 3.347, 5.287)Table 5
The optimal results of the QEPF-ATC approach.
Initial point
f
min
⁡
Iterations
Runtime (s)
(x
1, x
2, x
3, x
4, x
5, x
6, x
7)
1
2994.63432
42
272
(3.500, 0.699, 17, 7.300, 7.825, 3.392, 5.282)
2
2995.51327
38
215
(3.500, 0.703, 17, 7.305, 7.746, 3.353, 5.287)
3
2996.12356
30
198
(3.497, 0.696, 17, 7.300, 7.713, 3.326, 5.275)
4
2995.75363
17
97
(3.499, 0.700, 17, 7.299, 7.715, 3.349, 5.282)
5
2995.21554
52
306
(3.500, 0.700, 17, 7.300, 7.719, 3.350, 5.286)


## 5. Conclusion

Because ATC approach can flexibly organize the design optimization process according to the components of engineering systems, it has been widely applied in MDO. The QEPF-ATC approach is proposed to reduce the complicated coordination process and heavy computation cost faced by exiting ATC approaches, where quadratic exterior penalty function is adopted to be the replaced coordination strategy. Further, in order to deal with the design optimization of complex engineering systems, the QEPF-ATC and Kriging model combined approach is implemented, where Kriging models are used to replace simulation models with high accuracy in the solving process of MDO problems. The QEPF-ATC and Kriging model combined approach can avoid the impact of numerical noise brought by using simulation models during the iterative solving process and reduce computation cost. A geometric programming case is given to validate the usability and effectiveness of the QEPF-ATC approach. The QEPF-ATC and Kriging model combined approach has been successfully applied to solve the speed reducer design problem.It should be pointed out that using metamodels within the QEPF-ATC approach is not without its own set of potential drawbacks, which are well known and beyond the scope of this work. Prior to employing the metamodel assisted QEPF-ATC approach, the user should become familiar with the strengths and weaknesses associated with surrogate approximation. Variable-fidelity modeling provides a promising way to flexibly use approximate models in dealing with MDO problems. The combination of the QEPF-ATC approach and variable-fidelity modeling will be studied to design complex engineering systems in the future.



---

*Source: 685958-2015-03-25.xml*