# A Novel Hybrid Algorithm for Solving Multiobjective Optimization Problems with Engineering Applications

Authors: Lulu Fan; Tatsuo Yoshino; Tao Xu; Ye Lin; Huan Liu
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

An effective hybrid algorithm is proposed for solving multiobjective optimization engineering problems with inequality constraints. The weighted sum technique and BFGS quasi-Newton’s method are combined to determine a descent search direction for solving multiobjective optimization problems. To improve the computational efficiency and maintain rapid convergence, a cautious BFGS iterative format is utilized to approximate the Hessian matrices of the objective functions instead of evaluating them exactly. The effectiveness of the proposed algorithm is demonstrated through a comparison study, which is based on numerical examples. Meanwhile, we propose an effective multiobjective optimization strategy based on the algorithm in conjunction with the surrogate model method. This proposed strategy has been applied to the crashworthiness design of the primary energy absorption device’s crash box structure and front rail under low-speed frontal collision. The optimal results demonstrate that the proposed methodology is promising in solving multiobjective optimization problems in engineering practice.


---

## Body

## 1. Introduction

Multiobjective optimization problems (MOPs) are encountered in many fields, such as energy systems [1, 2], management [3], structural optimum design [4, 5], and scheduling [6, 7]. Generally, many problems can be modeled as multiobjective optimization problems, in which multiple conflicting objectives are to be optimized simultaneously. For example, in optimizing the design of a cycloid speed reducer, minimization of the volume and maximization of the efficiency should be considered simultaneously [8]. Vo-Duy et al. [9] attempted to minimize the weight of a laminated composite beam while maximizing its natural frequency. Optimization problems, which frequently appear in scientific research and engineering, are often MOPs, particularly in the design of structures [10, 11]. Yin et al. [12] studied the crashworthiness and reliability of a foam-filled bionic thin-walled structure based on bioinspired design. The optimization of a vehicle door structure with a hybrid material was investigated to achieve lightweight design [13].To obtain better performance, various multiobjective optimization strategies have been proposed and widely applied to engineering problems [14, 15]. Tire structure design was implemented based on multiobjective optimization [16]. A new multiswarm method was introduced, and its validity has been demonstrated on four benchmark functions and two structural problems [17]. A novel multiobjective particle swarm optimization algorithm was proposed based on decomposition to maintain diversified solutions [18]. However, most of these methods are evolutionary algorithms, which have the fatal drawback that their convergence cannot be theoretically guaranteed. In addition, the stopping criterion is difficult to determine according to a study by Coello [19]. Generally, a large number of iterations are preferred, which leads to considerable computation time and low efficiency, which cannot satisfy the engineering requirements [20].To solve the abovementioned problems, hybrid algorithms [21, 22] and a numerical algorithm based on Newton’s method [23] for MOPs were broadly proposed. For instance, for obtaining robust optimal solutions, a multiobjective optimization algorithm that combines particle swarm optimization and the gradient search technique was introduced [24], and Sun et al. [25] proposed a novel effective algorithm that aims to solve engineering design problems that involve uncertainties. To obtain a faster convergence rate and improve the stability of convergence and efficiency for solving engineering design problems, numerical algorithms gradually attracted increased attention from researchers [26, 27]. To handle vector problems, various techniques of multiobjective optimization have been studied, such as the ε-constraint method [28] and weighted max-min technique [29]. The scalarization technique, which was first introduced by Geoffrion [30], has been widely used to transform a multiobjective problem into a scalar one, and new insights into the weighted sum method have been reported in [31]. Yang et al. [32] have combined the weighted sum method and Newton’s method to develop a hybrid algorithm, which has been applied to practical engineering optimization.Generally, optimization design in a practical engineering problem is complicated and time-consuming. The surrogate model method has been properly utilized to simplify practical models in design problems [33, 34]. A surrogate-assisted evolution strategy was proposed and applied to multiobjective optimization [35]. Peng and Wang [36] established an effective adaptive surrogate for solving transfer trajectory optimization, the computing speed of which was almost 8 times faster than that of directly solving the actual model.In the past few years, thin-walled structures dramatically attracted the attention of researchers due to their excellent crashworthiness and extraordinarily light weight [37]. Thus, MOPs that were based on surrogate models that concerned thin-walled structures were investigated to improve the performance in terms of crashworthiness and light weight [38, 39]. Response surface methodology (RSM) was applied to the crashworthiness optimization of tapered thin-walled square tubes [40]. Sun et al. [41] have improved the energy absorption ability and realized the lightweight design of thin-wall structures based on multiobjective optimization and surrogate models.In this paper, a novel hybrid cautious BFGS quasi-Newton algorithm (denoted as CBQNA) is proposed based on the weighted sum technique and a cautious BFGS method for solving MOPs. To solve MOPs in practical engineering projects, an effective optimization strategy is proposed based on surrogated modeling and applied to crashworthiness optimization of the main absorption devices of the automobile. This paper has demonstrated that the proposed method is superior to evolutionary algorithm MOGA in obtaining uniformly distributed solutions and has faster calculation speed than MOGA and NSWFA.The rest of the paper is organized as follows. The hybrid algorithm CBQNA for MOP is proposed in Section2. The accuracy and efficiency of the novel algorithm are clearly demonstrated by comparison study on two benchmark functions in Section 3. Next, a rapid effective multiobjective optimization strategy and its engineering application are presented in Section 4. The study’s conclusions and future research are discussed in Section 5.


## 2. Fast Hybrid Algorithm for Multiobjective Optimization

In general, a multiobjective optimization problem with inequality constraints can be expressed as(1)min⁡Fx=F1x,F2x,…,FmxTs.t.gix≥0i=1,2,…,l,x∈D,where F(x)=(F1(x),F2(x),…,Fm(x))T is a vector function, m is the number of single-objective functions, and x=(x1,x2,…,xn)T is the decision variable vector, which consists of several variables that are generated in the feasible region D and satisfy the inequality constraint gi(x)≥0 for any i=1,2,…,l. To describe the optimal solutions of the MOP accurately and clearly, some definitions are clearly stated in the following section.Definition 1 (Pareto dominance).
Supposexv and xu are decision variables. xv is said to be dominated by xu (denoted as xv≺xu) if and only if, ∀j=1,2,…,m,Fj(xu) and Fj(xv) meet the condition(2)Fjxu<Fjxv∧Fjxu≤Fjxv.Definition 2 (Pareto-optimal solution).
Vectorxv is the optimal solution if and only if there is no solution x∈D that satisfies the dominance condition xv≺x.Definition 3 (Pareto-optimal set).
All Pareto-optimal solutions form the Pareto-optimal solution set (denoted asS) and it can be expressed as S={xv∈D∣x∈D∣,x≺xv}.Definition 4 (Pareto front).
The Pareto front (denoted asFS) consists of all corresponding values of function F(x) against all Pareto-optimal solutions xv∈S and(3)FS=Fxv=F1xv,F2xv,…,FmxvT∣xv∈S.


### 2.1. Treatment of Constraint Conditions by Penalty Function

To properly adapt the numerous effective unconstrained algorithms, the original problem (1) is converted into an unconstrained problem based on the Frisch interior penalty function. First, a penalty item is established as follows:(4)Px=-σk-1∑1lloggix,where 0<σk is a penalty coefficient, which is updated adaptively by (5) and gi(x) denotes the inequality constraint with i=1,2,…,l.(5)σk+1=ρkσk,where ρ>1 is a constant parameter. Then, the penalty function of the jth objective function Fj(xk) is constructed by(6)Qjx=Fjx+Px,where Fj(x) and P(x) are, respectively, defined in (1) and (4). Based on the above processes, problem (1) has been converted into an unconstrained problem, which can be described as(7)min⁡Fx=Q1x,Q2x,…,Qmxs.t.x∈D.To effectively solve the unconstrained problem, the preferred method is the BFGS quasi-Newton algorithm, which is superior to other similar methods, such as Newton’s method and the DFP method. The BFGS method is deduced from Newton’s method, which is based on Taylor expansion. The Taylor expansion of an unconstrained scalar function Qjx near xk can be expressed as(8)Qjxk+dk≈Qjxk+dkT∇Qjxk+12dkT∇2Qjxkdk,where ∇Qj(xk) and ∇2Qj(xk) are the gradient function and Hessian matrix, respectively, of target function Qj(xk) at xk. In quasi-Newton’s method, a symmetric matrix Bjk is utilized to approximate the Hessian matrix ∇2Qj(xk) instead of computing it directly. Next, inequality (8) can be converted to(9)Qjxk+dk≈Qjxk+dkT∇Qjxk+12dkTBjkdk.Consider ρk(d) is the function with respect to the variable d and the function ρk(d) is denoted as(10)ρkd=Qjxk+d-Qjxk=dT∇Qjxk+12dTBjkd.Assume that dk is the solution the minimum problem (11).(11)mind∈Rnρkd=Qjxk+d-Qjxk=dT∇Qjxk+12dTBjkd.If Bjk is a nonsingular matrix, the solution dk is obtained as (12) by solving problem (11).(12)dk=-Bjk-1∇Qjxk.The BFGS quasi-Newton’s method is an optimization method that uses (12) as its search direction. Generally, the positive-definite symmetric matrix Bjk is preferred to be similar in value to the Hessian matrix ∇2Qj(xk), which satisfies the equation in approximation (13).(13)∇Qj2xkSk≈Yjk,in which the vectors Yjk=∇Qjxk+1-∇Qjxk and Sk=xk+1-xk satisfy the condition(14)BjkSk=Yjk.Equation (14) is known as both the quasi-Newton equation and the quasi-Newton condition.


### 2.2. Weighted Sum Method

In this section, the weighted sum method is utilized to transform the vector functionF(x) into a scalar function, which is presented as(15)Gx=∑1mλjQjx,where λj∈(0,1) are weighting factors that meet the condition ∑j=1mλj=1 for j=1,2,…,m and λj can be divided into various sets according to the practical design demands. Based on the abovementioned processes, the original multiobjective problem with inequality constraints has been converted to an unconstrained scalar problem, which is stated as(16)min⁡Gxx∈D.


### 2.3. Hybrid Algorithm for MOP

In this section, a cautious BFGS quasi-Newton algorithm (CBQNA) for the MOP based weighted sum method is proposed. Consider the vector function that is shown in problem (1), which has been translated into a single-objective optimization problem (16). First, based on the interior penalty function method and the weighted sum technique, a quasi-Newton direction dk for the MOP is obtained by solving(17)dk=-∑j=1mλjBjk-1∑jmλj∇Qjxk.In this paper, the approximation matrix Bjk is updated by a cautious BGFS equation (18), which was introduced in [42] and can guarantee the positive definiteness of Bjk after each iteration.(18)Bjk+1=Bjk-BjkSkSkTSkTBjkSk+YjkYjkTSkTYjkifYjkTSkSk2≤ε∇QjxkγBjkelse,where ε and γ are positive constant parameters; the settings in [43] are followed to maintain high accuracy. The search direction dk that is obtained by solving (17) is always a descent direction when the condition SkTYjk≥0 is satisfied. If this condition is not satisfied, the gradient direction will be replaced by another search direction, which is expressed as(19)dk=-∑jmλj∇Qjxk.Then, the classical Armijo rule is used to determine a proper step length αk along quasi-Newton’s search direction dk. In general, the Armijo rule for search direction dk of function f(x) is(20)fx+αkdk≤fx+βαk∇fxTdk,where β∈(0,1) is a parameter. Under the condition of maintaining the descent search direction, a full step length is preferred. However, the value αk=1 often does not satisfy inequality (20), so the largest αk that satisfies this inequality is always selected for general cases. Furthermore, the positive definiteness of Bk can be preserved by utilizing updating equation (18) and setting β∈(0,1/2), and the full step length along the descent search direction is usually obtained according to the study [44]. Thus, in this paper, we rewrite the Armijo condition as(21)Gx+αkdk≤Gx+μαk∑1mλj∇QjxTdk,where G(x) is expressed in (15) and μ∈(0,1) is a constant parameter.Therefore, based on the previous sections, the hybrid algorithm CBQNA can be implemented by following Step  0 to Step  7.Step  0 (initialization). Choose an initial value x0∈D⊂Rn, a symmetric positive-definite matrix B0∈Rn×n, and a sufficiently small constant ε1>0. Moreover ρ>1,0<σ0<1,μ∈(0,1), and γ,ε>0.Step  1 (stopping criterion 1). If xk+1-xk∨dk≤ε1, then stop; else, go to Step  2.Step  2 (descent search direction). The search direction dk is obtained by solving (17) or (19), according to whether the descent condition is satisfied or not.Step  3 (step length). Compute a step length αk along the search direction dk such that the Armijo-rule inexact line search rules (21) are satisfied.Step  4 (iterative procedure). Set xk+1=xk+αkdk. If xk+αkdk∈D,k:=k+1; otherwise, compute xk+1 by(22)xk+1=xLifxk+αkdk<xLxUelsexk+αkdk>xU.Step  5 (stopping criterion 2). If P(xk)≤ε1 stop with x∗=xk; else, set σk=ρkσk and return to Step  1.


### 2.4. Analysis of the Proposed Algorithm

The novel algorithm that is proposed in this study is sensitive to the initial value, similar to numerical methods that are based on gradient information. However, this disadvantage can be overcome by repeatedly adjusting the starting point or selecting some preferred value that is obtained by other methods. To achieve globally optimal solutions by CBQNA, the Hessian matrices of the optimization functions are approximated by a cautious BFGS algorithm in this study. The proposed algorithm is primarily utilized to solve MOPs in engineering projects, the mathematical models of which are usually replaced by approximate convex functions that can ensure that the global optimal solution of the problem be achieved. This method has a fast convergence rate and can obtain optimal solutions in few iterations, which is always difficult for evolutionary algorithms. Particularly, the proposed algorithm is more effective in solving engineering problems than algorithms such as Newton’s method, because it avoids computing second-order information of the objective functions, which does not exist in the general engineering case.


## 2.1. Treatment of Constraint Conditions by Penalty Function

To properly adapt the numerous effective unconstrained algorithms, the original problem (1) is converted into an unconstrained problem based on the Frisch interior penalty function. First, a penalty item is established as follows:(4)Px=-σk-1∑1lloggix,where 0<σk is a penalty coefficient, which is updated adaptively by (5) and gi(x) denotes the inequality constraint with i=1,2,…,l.(5)σk+1=ρkσk,where ρ>1 is a constant parameter. Then, the penalty function of the jth objective function Fj(xk) is constructed by(6)Qjx=Fjx+Px,where Fj(x) and P(x) are, respectively, defined in (1) and (4). Based on the above processes, problem (1) has been converted into an unconstrained problem, which can be described as(7)min⁡Fx=Q1x,Q2x,…,Qmxs.t.x∈D.To effectively solve the unconstrained problem, the preferred method is the BFGS quasi-Newton algorithm, which is superior to other similar methods, such as Newton’s method and the DFP method. The BFGS method is deduced from Newton’s method, which is based on Taylor expansion. The Taylor expansion of an unconstrained scalar function Qjx near xk can be expressed as(8)Qjxk+dk≈Qjxk+dkT∇Qjxk+12dkT∇2Qjxkdk,where ∇Qj(xk) and ∇2Qj(xk) are the gradient function and Hessian matrix, respectively, of target function Qj(xk) at xk. In quasi-Newton’s method, a symmetric matrix Bjk is utilized to approximate the Hessian matrix ∇2Qj(xk) instead of computing it directly. Next, inequality (8) can be converted to(9)Qjxk+dk≈Qjxk+dkT∇Qjxk+12dkTBjkdk.Consider ρk(d) is the function with respect to the variable d and the function ρk(d) is denoted as(10)ρkd=Qjxk+d-Qjxk=dT∇Qjxk+12dTBjkd.Assume that dk is the solution the minimum problem (11).(11)mind∈Rnρkd=Qjxk+d-Qjxk=dT∇Qjxk+12dTBjkd.If Bjk is a nonsingular matrix, the solution dk is obtained as (12) by solving problem (11).(12)dk=-Bjk-1∇Qjxk.The BFGS quasi-Newton’s method is an optimization method that uses (12) as its search direction. Generally, the positive-definite symmetric matrix Bjk is preferred to be similar in value to the Hessian matrix ∇2Qj(xk), which satisfies the equation in approximation (13).(13)∇Qj2xkSk≈Yjk,in which the vectors Yjk=∇Qjxk+1-∇Qjxk and Sk=xk+1-xk satisfy the condition(14)BjkSk=Yjk.Equation (14) is known as both the quasi-Newton equation and the quasi-Newton condition.


## 2.2. Weighted Sum Method

In this section, the weighted sum method is utilized to transform the vector functionF(x) into a scalar function, which is presented as(15)Gx=∑1mλjQjx,where λj∈(0,1) are weighting factors that meet the condition ∑j=1mλj=1 for j=1,2,…,m and λj can be divided into various sets according to the practical design demands. Based on the abovementioned processes, the original multiobjective problem with inequality constraints has been converted to an unconstrained scalar problem, which is stated as(16)min⁡Gxx∈D.


## 2.3. Hybrid Algorithm for MOP

In this section, a cautious BFGS quasi-Newton algorithm (CBQNA) for the MOP based weighted sum method is proposed. Consider the vector function that is shown in problem (1), which has been translated into a single-objective optimization problem (16). First, based on the interior penalty function method and the weighted sum technique, a quasi-Newton direction dk for the MOP is obtained by solving(17)dk=-∑j=1mλjBjk-1∑jmλj∇Qjxk.In this paper, the approximation matrix Bjk is updated by a cautious BGFS equation (18), which was introduced in [42] and can guarantee the positive definiteness of Bjk after each iteration.(18)Bjk+1=Bjk-BjkSkSkTSkTBjkSk+YjkYjkTSkTYjkifYjkTSkSk2≤ε∇QjxkγBjkelse,where ε and γ are positive constant parameters; the settings in [43] are followed to maintain high accuracy. The search direction dk that is obtained by solving (17) is always a descent direction when the condition SkTYjk≥0 is satisfied. If this condition is not satisfied, the gradient direction will be replaced by another search direction, which is expressed as(19)dk=-∑jmλj∇Qjxk.Then, the classical Armijo rule is used to determine a proper step length αk along quasi-Newton’s search direction dk. In general, the Armijo rule for search direction dk of function f(x) is(20)fx+αkdk≤fx+βαk∇fxTdk,where β∈(0,1) is a parameter. Under the condition of maintaining the descent search direction, a full step length is preferred. However, the value αk=1 often does not satisfy inequality (20), so the largest αk that satisfies this inequality is always selected for general cases. Furthermore, the positive definiteness of Bk can be preserved by utilizing updating equation (18) and setting β∈(0,1/2), and the full step length along the descent search direction is usually obtained according to the study [44]. Thus, in this paper, we rewrite the Armijo condition as(21)Gx+αkdk≤Gx+μαk∑1mλj∇QjxTdk,where G(x) is expressed in (15) and μ∈(0,1) is a constant parameter.Therefore, based on the previous sections, the hybrid algorithm CBQNA can be implemented by following Step  0 to Step  7.Step  0 (initialization). Choose an initial value x0∈D⊂Rn, a symmetric positive-definite matrix B0∈Rn×n, and a sufficiently small constant ε1>0. Moreover ρ>1,0<σ0<1,μ∈(0,1), and γ,ε>0.Step  1 (stopping criterion 1). If xk+1-xk∨dk≤ε1, then stop; else, go to Step  2.Step  2 (descent search direction). The search direction dk is obtained by solving (17) or (19), according to whether the descent condition is satisfied or not.Step  3 (step length). Compute a step length αk along the search direction dk such that the Armijo-rule inexact line search rules (21) are satisfied.Step  4 (iterative procedure). Set xk+1=xk+αkdk. If xk+αkdk∈D,k:=k+1; otherwise, compute xk+1 by(22)xk+1=xLifxk+αkdk<xLxUelsexk+αkdk>xU.Step  5 (stopping criterion 2). If P(xk)≤ε1 stop with x∗=xk; else, set σk=ρkσk and return to Step  1.


## 2.4. Analysis of the Proposed Algorithm

The novel algorithm that is proposed in this study is sensitive to the initial value, similar to numerical methods that are based on gradient information. However, this disadvantage can be overcome by repeatedly adjusting the starting point or selecting some preferred value that is obtained by other methods. To achieve globally optimal solutions by CBQNA, the Hessian matrices of the optimization functions are approximated by a cautious BFGS algorithm in this study. The proposed algorithm is primarily utilized to solve MOPs in engineering projects, the mathematical models of which are usually replaced by approximate convex functions that can ensure that the global optimal solution of the problem be achieved. This method has a fast convergence rate and can obtain optimal solutions in few iterations, which is always difficult for evolutionary algorithms. Particularly, the proposed algorithm is more effective in solving engineering problems than algorithms such as Newton’s method, because it avoids computing second-order information of the objective functions, which does not exist in the general engineering case.


## 3. Numerical Examples

In this part, the performance of the proposed hybrid algorithm is compared with those of the Multiobjective Genetic Algorithm (MOGA) and a hybrid algorithm (denoted as NSWFA), which was proposed in [32], on two benchmark tests. To evaluate this novel algorithm, the CPU time and the number of iterations for obtaining a Pareto-optimal solution are recorded. All MATLAB procedures are run in an identical computing environment. Each test problem was solved 100 times with initial points that were generated from a uniform random distribution over the feasible zone. All attempts to solve the two examples were terminated when the stopping criterion was satisfied or 2000 iterations were performed.Example 1.
A simple biobjective problem with one variable, which was derived from a MATLAB tool box [45]:(23)minF1x=x2F2x=x-22s.t.-1000≤x≤1000.Example 2.
A biobjective problem with two variables, which was obtained from [46]: (24)minF1x=x14+x24+x1x2-10x14-x12x22F2x=x14+x24+x1x2-x12x22s.t.-5≤x1≤5;-5≤x2≤5.According to the research in [32], the Pareto-optimal solutions that are obtained by algorithm NSWFA are more uniform than those by MOGA. For simplicity, only the Pareto-optimal fronts for the two tests that are obtained by NSWFA and CBQNA are shown in Figures 1 and 2, respectively. In addition, the characteristics of the three algorithms in terms of the number of iterations, CPU time, and convergence condition are intuitively presented in Table 1.Table 1
Detailed comparison of the three algorithms in terms of computational efficiency.
Algorithm
Stopping criterion
CPU time
Iterations
Convergence
Test 1
MOGA
10−6
391.9 s
2000
No
NSWFA
10−6
5.250 s
615
Yes
CBQNA
10−6
0.423 s
200
Yes
Test 2
MOGA
10−5
391.9 s
2000
No
NSWFA
10−5
13.23 s
730
Yes
CBQNA
10−5
2.600 s
715
YesFigure 1
Pareto-optimal fronts of Test 1, obtained by CBQNA and NSWFA.Figure 2
Pareto-optimal fronts of Test 1, obtained by CBQNA and NSWFA.Figures1 and 2 illustrate that the proposed algorithm performs similarly to NSWFA in terms of distribution uniformity of the solutions on Test 1, while the performance of CBQNA is slightly worse than that of NSWFA on Test 2. The Pareto-optimal fronts that are obtained by NSWFA and the proposed algorithm CBQNA are closely related to the initial value. This finding is observed because the two algorithms behave similarly when searching for the optimal solutions. In addition, according to Table 1, CBQNA is superior to the other two algorithms on the two examples, in either the number of iterations or CPU time. In addition, when the solution accuracy is the same among all three algorithms, the proposed algorithm uses the minimal number of iterations and has the fastest convergence rate. In summary, the proposed CBQNA performs nearly as well as NSWFA and outperforms MOGA in terms of obtaining evenly distributed Pareto-optimal solutions, which may be worse than the results obtained through other methods. Furthermore, the computational efficiency of the novel algorithm in this study highlights the advantages of CBQNA over MOGA and NSWFA, which are clearly shown in Table 1. Particularly, when determining the descent search direction, the proposed algorithm avoids computing the Hessian matrices of objective functions, which is of great significance for solving practical engineering problems.


## 4. Effective Optimization Strategy and Its Practical Engineering Application

### 4.1. Optimization Strategy for MOPs

To reduce the complexity of engineering multiobjective optimization problems, a valid optimization strategy that is based on surrogate modeling and the proposed algorithm CBQNA is properly introduced. The effective strategy was established according to the following four processes.


#### 4.1.1. Establishment of the Mathematical Model

To effectively solve optimization design problems, the optimization objectives and design variables need to be determined correctly. The preferred targets are often related to properties of the optimized structures. Analysis of sensitivity is usually applied when choosing decision variables, which greatly affects the targets. The design domains and constraint conditions are always determined according to the limiting dimension and strength requirements of the structure, scheme, or system. Hence, an accurate mathematical expression for the multiobjective optimization problem should be constructed based on the aforementioned analysis to solve the MOPs.


#### 4.1.2. Design of Experiments

The establishment of an agent model with high precision is based on a good design scheme. The setting of reasonable factor levels and variable span play a significant role in the experimental design, especially in multifactor and multiobjective problems. Common methods for experimental design are orthogonal design and central composite design (CCD). In addition, Taguchi robust design has been applied in experimental design for its better performance [25]. To describe how to construct sample points and determine the optimal number of samples by orthogonal design and CCD, the procedural details of the two design methods can be separately summarized as follows.(1) Orthogonal Design. In general, the format of the orthogonal design is comparatively fixed if the number and the levels of the design variables have been given. That is, the optimal number of samples and its corresponding experimental factors can be determined if the ranges of the design and optimal levels are provided in advance. The distinctive features of schemes that are planned by orthogonal design are orthogonality and even distribution. Furthermore, orthogonal design requires fewer experiments in establishing a surrogate model, which improves the efficiency.(2) Central Composite Design (CCD). A complete scheme that is obtained by central composite design has three different experimental patterns. For example, if there are three design variables, the distribution of the sample points that are selected by CCD is illustrated in Figure 3. The eight points that are located at the summit of the cube are of the same experimental format, the six sample points that are located at the end of the center axis are of the same type, and the center point in red is a sample that needs to be repeated six times. Thus, for three variables, the optimal number of samples is twenty.Figure 3
Distribution of samples planned by CCD.


#### 4.1.3. Establishment of Surrogate Model

Surrogate modeling has been broadly applied to many fields for its simple theory and excellent behavior. In the general case, the response surface methodology (RSM) is preferred by researchers for its simplicity and validity. Particularly, the RSM is more effective when dealing with multiobjective problems with multiple variables. Therefore, considering the characteristics of the proposed algorithm, the establishment of a linear, quadratic, or cubic polynomial surface is favored for describing the relationships between objective responses and design variables.Generally, a linear polynomial model is defined as follows:(25)y^=a+b1x1+b2x2+⋯+bmxm,where y^ is an approximate expression of objective or constraint functions, which is known as the response surface; xi is design variable i, with i=1,2,…,m; and m is the number of variables. The constant coefficients a,b1,b2,…,bm in (25) can be calculated in terms of the principle of least squares. Similarly, a quadratic polynomial model is usually constructed as in expression (26) and the cubic polynomial surface can be obtained in the same manner.(26)y^=a+∑1mbixi+∑1mbiixi2+∑i<kbikxixk.For simplicity, expression (26) can be converted to a linear form as follows:(27)y^=a+∑i=1WciXi,where W=m(m+3)/2 and X=(X1,X2,…,Xm) denote items with x1,x2,…,xm,x12,x22,…,xm2 and x1x2,x1x3,…,xm-1xm, respectively.To accurately estimate the precision of the approximated models, the correlation coefficientR2 and adjusted index Radj2 need to be calculated according to formulas (28) and (29), respectively.(28)R2=∑i=1My^i-y¯i2∑i=1Myi-y¯i2(29)Radj2=1-∑i=1Myi-y^i2M-1∑i=1Myi-y¯i2M-k-1,where M denotes the number of factors and y^i, y¯i  are the predictive value and the actual value, respectively, of the average measurement of y. The closer R2 and Radj2 are to 1, the better the surrogate model approximates the original problem. If they do not meet the accuracy requirement, the experimental points need to be reset to establish another model until the precision is satisfied.


#### 4.1.4. Pareto Solutions Obtained by CBQNA

In this process, the approximate mathematical expressions of the primal multiobjective problems are first translated into a standard normalized multiobjective problem. Next, the Pareto-optimal solution set of the normalized problem is directly obtained by the CBQNA algorithm following the steps that are described in Section2. In practical projects, each Pareto solution corresponds to a design scheme, which can provide references for designers to optimize the original design.Based on these main steps, an effective optimization strategy for solving multiobjective optimization problems in engineering projects is in Figure4.Figure 4
Effective optimization strategy based on CBQNA.


### 4.2. Crashworthiness Optimization of Crash Box

With the boom of the automobile industry, the importance of passive safety has attracted increasing attention from researchers. To improve the energy absorption abilities of the main absorbers, such as the crash box, the longitudinal beam and door structure are vital to the protection of occupants. Thus, the analysis of the crashworthiness of the automobile energy absorbing box is essential and crucial.In this section, the proposed optimization strategy is applied to crashworthiness optimization of an automobile crash box. Consider a crash box that is composed of four thin-walled plates, as shown in Figure5. This thin-walled structure has been optimized by Yang et al. in study [32], but the impact condition that was used cannot correctly simulate the actual operation. Thus, a better simulation model under frontal collision is built in this paper.Figure 5
Bumper-crash box structure of the automobile.To effectively simulate the performance of the crash box under a low-speed impact, one end of each crash box has been fully constrained, while the other end is subject to front collision by a rigid wall of infinite mass at a low speed of 4 m/s along thex-axis, as shown in Figure 5. The original design of the crash is shown in Figure 6. The total impact time is 0.05 s.Figure 6
Original design of the crash box.


#### 4.2.1. Mathematical Model for the MOP of the Crash Box Structure

The energy absorption and the maximal collision force, as the two major indices of the crashworthiness ability of the crash box, have been considered simultaneously. The crash box is a thin-walled structure, whose properties are affected dramatically by its wall thickness. In addition, when designing thin-walled structures, especially for automobile absorbers, the maximal crushing force and the mass should each be constrained below a maximum value. Thus, improving energy absorptionEin and reducing maximal impact force Fmax have been selected as simultaneous objectives, while the thicknesses of the four plates are chosen as design variables. Moreover, the mass of the entire structure and the maximal impact force are constrained by upper limits. Consequently, an MOP of the crash box that is based on crashworthiness is mathematically expressed as(30)min⁡F1=1Einx,F2=Fmaxxs.t.Fim≤F0;Mt≤M0⌻x∈xL,xU,where F1(x) and F2(x) are two objective functions of the vector x=(t1,t2,t3,t4), which consists of four decision variables, which are subject to boundary condition ti∈[1,3]mm,i=1,2,3,4. The conditions Fmax≤F0 and Mt≤M0 constrain the impact force and the total mass should be maintained below an upper value. In this paper, we set F0=174 kN and M0=4.2 kg.


#### 4.2.2. Experiment Design

The thicknesses of the four thin-walled plates are concurrently considered by DOE. The combined orthogonal experimental method is utilized to construct a 17-group experimental dataset in the design domain.


#### 4.2.3. Construction of Response Surface Models

The response surface methodology is utilized to construct surrogate models instead of the original complicated problem. The mathematical expressions for the two objectives in terms of the four variables are stated as (31) and (32). In addition, the approximate mode of the constraint function is presented in expression (33).(31)Einx=-3397.09+2967.8t1+3368.16t2+1300.84t3+254.58t4-157.9t12-211.4t22-31.41t32-38.69t42-684.18t1t2-93.37t1t3-307.1t2t3+85.12t2t4-65.4t3t4-38.7t1t4;(32)Fmaxx=-76.71-13.62t1-10.5t2+107.4t3+142.8t4+28.14t12+24.44t22-18.43t32-18.38t42-11.13t1t2-25.40t1t3-6.38t2t3-39.58t2t4+13.99t3t4-10.96t1t4;(33)Mtx=-0.0128+0.389t1+1.19t2+0.346t3+0.183t4.After calculation, the coefficients of determination of energy absorption, maximal impact force, and mass constraint are 98.60%, 99.11%, and 1, respectively, while the adjusted coefficients of determination are 95.6%, 97.11%, and 1, respectively. The calculated indicators demonstrate that the approximate models of the two objective and constraint functions have excellent precision.


#### 4.2.4. MOP by CBQNA

The hybrid algorithm CBQNA, which is introduced in Section2, is applied to solve the optimization design problem. For adapting the surrogate model to the novel algorithm, a normalization process is implemented first. Thus, problem (30) can be converted into a standard mathematical problem, which is expressed as follows:(34)minF1x=EmaxEinxF2x=FmaxxFas.t.1mm≤ti≤3mm,i=1,2,3,4.Fmaxx≤174kN;Mtx≤4.2kg,where Emax and Fa are the maximal value of the energy absorption and average value of the maximal collision force, respectively.After optimization, 100 groups of solutions of the optimization problem are obtained, which takes only 5.6 seconds on average. The Pareto-optimal front, which is composed of 100 Pareto points, is presented in Figure7, and the relationship between energy absorption and maximal impact force is depicted in Figure 8.Figure 7
Pareto-optimal front obtained by CBQNA.Figure 8
Relationship between the energy absorption and the maximal impact force.In terms of the properties of energy absorption ability and maximal collision force, the optimal design is specified byt1=2.0,t2=1.98,t3=1.98, and t4=1.98. Regarding the lightweight design, the mass of the thin-walled plates is also a significant factor in selecting the best design. The relationship between the energy absorption, maximal impact force, and total mass of the bumper absorbing box is depicted in Figure 9. The comparison between the original design and the optimal result in terms of energy absorption, maximal collision force, and mass is presented in Table 2.Table 2
Comparison of design variables and the performances before and after optimization.
t
i
/
m
m
E
i
n
/
J
F
m
a
x
/
k
N
M
a
/
k
g
t
1
t
2
t
3
t
4
Before
1.4
1.8
2.3
1.8
5605.2
173.25
3.83
After
1.86
1.78
2.0
1.96
5777.8
146.64
3.89Figure 9
Objective responses and total mass of the crash box structure.Table2 shows that, compared with the original design, the performance improves greatly after optimization. The energy absorption improves by 3.08% and the maximal impact force reduces by 14.87%, while the total mass increases slightly by 1.57%. The collision conditions of the crash box at 0.01 s, 0.02 s, 0.03 s, and 0.04 s are presented in Figure 10. From the collision deformation, the optimized crash box structure achieves an ideal deformation and satisfies the practical engineering requirements.Figure 10
Deformation of the crash box at different times.


### 4.3. Optimization Design of an S-Shaped Thin-Walled Beam

The frontal rail that is connected with the crash box is also a main deformation and energy absorption structure in a frontal crash, which is shown in Figure11.Figure 11
Simplified model of an S-shaped front rail structure.Generally, the vehicle front rail is a thin-walled structure whose deformation and energy absorption characteristics are similar to those of a general rectangular thin-walled beam. Thus, the crashworthiness of the thin-walled beam has dramatic significance to the frontal crash safety of the automobile. In this paper, an attempt has been made to study an S-shaped simplified model under frontal impact conditions.


#### 4.3.1. Establishment of the Finite-Element Model

Firstly, the geometry model of S-shaped structure is established with Catia, then its finite-element model is simulated by Hypermesh, and the simulation process is carried by LS-DYNA. The FE model of the S-shaped thin-walled beam with rectangular cross-section is displayed in Figure12.The S-shaped beam is impacted by a rigid wall that weights 450 kg at the speed of 10 m/s. As shown in Figure 13, the rear end of the beam is fully constrained and the other end is allowed to move only along the X-direction. The material that is utilized for the S-shaped beam is mild steel with the following material properties: Young’s modulus = 206 Gpa, Poisson’s ratio = 0.3, and density = 7800 kg/m3. The relationship between stress and strain rate (shown in Figure 14) is also considered.Figure 12
FE model of the S-shaped beam.Figure 13
Boundary condition for S-shaped front rail.Figure 14
Typical engineering stress-strain curve of the material [47].


#### 4.3.2. Mathematical Model for MOP of the Front Rail

The cross-sectional dimensions of the thin-walled structure are the main factor that affects the properties of energy absorption and maximal impact resistance. The widthw, height h, and thickness t of the cross-section have been chosen as design variables in this study, and the cross-section of an S-shaped thin-walled structure is shown in Figure 15.Figure 15
Cross-section of an S-shaped beam.In addition, energy absorption and maximal impact force are major evaluation indices of the crashworthiness ability of the S-shaped front rail beam.Hence, the mathematical model for multiobjective optimization of the front rail, which consists of two objective functions in terms of three design variables, is constructed as(35)min⁡Fx=F1x=1ExF2x=Fxs.t.hmin≤h≤hmaxwmin≤w≤wmaxtmin≤t≤tmax,where x=(w,h,t)T is the vector of design variables; F(x) is a vector function that consists of two objectives, namely, F1(x) and F2(x); hmin and hmax denote the lower and upper limits of the height of the cross-section; wmin and wmax denote the lower and upper limits of the width; and tmin and tmax are the lower and upper limits of the thickness of the wall, respectively. In this study, to satisfy the limits of the structure dimensions and the strength demands, hmin=136mm,hmax=146mm,wmin=56mm,wmax=68mm,tmin=1.4mm, and tmax=2.2mm.


#### 4.3.3. Design of Experiments

Following the novel optimization strategy’s procedures, orthogonal design is utilized to generate 25 design sampling points (i.e., four levels for each of the three design variables), which are clearly presented in Table3 with their corresponding responses.Table 3
Data for 25 groups, obtained by orthogonal experimental design.
Number
Factors
E  (J)
F  (kN)
w
h
t
(1)
58
136
1.4
11444.7
182.067
(2)
58
139
1.6
12943.7
214.618
(3)
58
142
1.8
14335.4
235.309
(4)
58
145
2.0
15644
251.437
(5)
58
148
2.2
16845.9
287.804
(6)
60
136
1.6
12910.2
212.551
(7)
60
139
1.8
14301
233.593
(8)
60
142
2.0
15603.6
248.585
(9)
60
145
2.2
16805.2
272.659
(10)
60
148
1.4
12257.2
196.856
(11)
62
136
1.8
14242.8
221.679
(12)
62
139
2.0
15547.9
246.477
(13)
62
142
2.2
16643.2
299.146
(14)
62
145
1.4
12236.6
197.617
(15)
62
148
1.6
13765.4
231.555
(16)
64
136
2.0
15383.9
259.826
(17)
64
139
2.2
16578.2
299.329
(18)
64
142
1.4
12200.4
194.745
(19)
64
145
1.6
13730.5
229.025
(20)
64
148
1.8
15148.8
241.987
(21)
66
136
2.2
16513.7
299.543
(22)
66
139
1.4
12157.1
195.991
(23)
66
142
1.6
13679.3
224.678
(24)
66
145
1.8
15094.3
234.490
(25)
66
148
2.0
16400.2
272.619


#### 4.3.4. Agent Model and Precision Analysis

A surrogate model of two objectives for crashworthiness optimization for the frontal rail is also established to lower the complexity of the original problem. According to the 25 plans of the DOE, which are shown in Table3, the quadratic response surface models of the energy absorption and maximal impact force responses were constructed in expressions (36) and (37), respectively.(36)E′x=-70973.4+958.656w+614.173h-3.827wh+57.023wt+63.274ht-3.582w2-1.4842h2-1890.174t2(37)F′x=111.854-20.210w+7.560h-0.025wh+2.731wt-1.479ht+0.165w2-0.009h2+43.110t2.The surrogate model’s fitness was evaluated based on coefficientR2 and its adjusted form Radj2. The values of coefficient R2 for the energy absorption and collision force reach 99% and 98%, respectively, and the correction coefficient Radj2 values are 99% and 95%, respectively. The result indicates that the approximation ability of this surrogate model is very good.


#### 4.3.5. Optimization Design

To improve the efficiency of optimization, the hybrid algorithm CBQNA is utilized to minimize the two response functions simultaneously. First, based on the new algorithm, a normalization procedure is applied to unify the magnitudes of the two responses functions. Thus, problem (35) is converted into a standard biobjective optimization problem with inequality constraints, which is expressed as follows:(38)min⁡Fx=F1x=EaE′xF2x=F′xFa(39)s.t.h-136146-h≥0;w-5866-w≥0,t-1.42.2-t≥0;where Ea=14496.5J and Fa=234.9 kN are the average value of energy absorption and the maximal collision force for the 25 groups, respectively.


#### 4.3.6. Optimization Results Solved by CBQNA

Problem (38) is solved by the new proposed algorithm with random initial value in accordance with constraints. To maintain the diversity of design selection, 200 groups of optimal solutions are acquired regularly by adjusting weighting factors. The Pareto-optimal front, which is constructed with 200 solutions, is depicted in Figure 16. In addition, the Pareto-optimal front of the front rail collision problem is shown in Figure 17. When determining the dimensions of the front rail, we should focus on not only the crashworthiness ability but also the weight requirement of the automobile. Thus, the mass of the thin-walled beam is also considered in evaluating the crashworthiness performance in Table 4.Table 4
Optimal solutions and their corresponding objective responses.
Design variables
E
(
J
)
F
(
k
N
)
M
(
k
g
)
w
h
t
(1)
63.18507
145.4626
2.1208
17081.0
272.889
7.19
(2)
63.23390
145.5450
1.8717
15538.6
245.962
6.36
(3)
63.23933
145.5548
2.0046
16087.7
263.545
6.80
(4)
63.23917
145.5540
1.9672
16087.7
263.545
6.67
(5)
63.17605
145.4450
2.1677
17048.5
267.451
7.36
(6)
63.05213
145.2380
1.8668
15485.6
245.663
6.32
(7)
63.10070
145.3174
1.7284
14311.9
233.394
5.89
(8)
63.16113
145.4234
1.8880
15506.7
245.863
6.43
(9)
63.22119
145.5222
1.9085
15511.4
245.155
6.46
(10)
63.21939
145.5207
1.8563
15509.1
243.437
6.33
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
(71)
63.20308
145.4854
1.8738
15506.2
245.319
6.38
(72)
63.14951
145.3912
1.7064
14323.5
233.614
5.81
(73)
63.18035
145.4439
2.1700
17036.8
276.266
7.38
(74)
63.07671
145.2661
1.9681
16030.8
257.892
6.68
(75)
63.18094
145.4458
1.9650
16049.1
258.247
6.69
(76)
63.14457
145.3848
2.0132
16037.5
260.135
7.02
(77)
63.19502
145.4693
2.1381
16586.5
266.359
7.26
(78)
63.20000
145.4776
1.8571
15677.5
240.558
6.23
(79)
63.11433
145.3301
2.1676
17037.7
265.008
7.35
(80)
63.20090
145.4790
1.9712
16056.5
257.411
6.70
(81)
63.19128
145.4618
2.0219
16054.8
264.999
6.86
(82)
63.15801
145.4045
1.9299
15504.9
246.552
6.56
(83)
63.12779
145.3543
1.9100
15507.3
245.189
6.49
(84)
63.20623
145.4870
1.9823
16052.2
258.096
6.73
(85)
63.18679
145.4535
1.9746
16055.6
257.415
6.70
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…Figure 16
Pareto-optimal front obtained by CBQNA.Figure 17
Relationship between the energy absorption and the maximal impact force.According to the properties of the best plans, the cross-sectional dimensions are selected, and the thickness of the wall is set tox=63.2,146,1.8mm, as shown in line 78 in Table 4. After the simulation by software LS-DYNA, the energy absorption and collision force responses, the specific energy absorption (SEA), and the mass of the structure, compared with the original design, are presented in Table 5. Figure 18 illustrates the deformation of the S-shaped frontal rail thin-walled structure changes over time.Table 5
Design variables and response values before and after optimization.
Design variables (mm)
E
(
J
)
F
(
k
N
)
M
(
k
g
)
w
h
t
Original Design
62.0
142
2.00
15616.8
247.9
6.45
Optimized Result
63.2
146
1.80
15077.5
240.5
6.03Figure 18
Deformation of the S-shaped beam changes over time.Although, compared with the original design, the internal energy absorption of the optimal design is decreased slightly, the maximal impact and the mass are decreased by 3.0% and 6.5%, respectively. Moreover, the SEA, which is the main indicator of the energy absorption ability, drops by 3%, which signifies that the crashworthiness of the frontal rail is improved. Therefore, all performance indicators of the optimal result meet the requirements.


### 4.4. Analysis of Optimal Results

The optimization strategy that is based on CBQNA has been applied to the optimization design of the crash box and the S-shaped front rail. The optimized results demonstrate that the optimization strategy that is based on the CBQNA algorithm is valuable in solving MOPs in engineering projects. To evaluate the efficiency of the proposed method, the two engineering examples in this paper are each executed 30 times. The calculation results show that it takes approximately 5.6 s to optimize the crash box by this novel algorithm, compared to approximately 198 s by NSWFA [32]. Furthermore, the constraint condition of the crash box that is simulated in this paper is more consistent with the actual frontal collision than that in [32]. The calculation time of the optimization of the S-shaped front rail is approximately 32.8 s for 200 groups’ optimal solutions, which demonstrates that the optimization strategy that is based on the proposed algorithm is effective and promising in solving complicated engineering problems in which the objective function’s second derivative is difficult to calculate or does not exist.


## 4.1. Optimization Strategy for MOPs

To reduce the complexity of engineering multiobjective optimization problems, a valid optimization strategy that is based on surrogate modeling and the proposed algorithm CBQNA is properly introduced. The effective strategy was established according to the following four processes.


### 4.1.1. Establishment of the Mathematical Model

To effectively solve optimization design problems, the optimization objectives and design variables need to be determined correctly. The preferred targets are often related to properties of the optimized structures. Analysis of sensitivity is usually applied when choosing decision variables, which greatly affects the targets. The design domains and constraint conditions are always determined according to the limiting dimension and strength requirements of the structure, scheme, or system. Hence, an accurate mathematical expression for the multiobjective optimization problem should be constructed based on the aforementioned analysis to solve the MOPs.


### 4.1.2. Design of Experiments

The establishment of an agent model with high precision is based on a good design scheme. The setting of reasonable factor levels and variable span play a significant role in the experimental design, especially in multifactor and multiobjective problems. Common methods for experimental design are orthogonal design and central composite design (CCD). In addition, Taguchi robust design has been applied in experimental design for its better performance [25]. To describe how to construct sample points and determine the optimal number of samples by orthogonal design and CCD, the procedural details of the two design methods can be separately summarized as follows.(1) Orthogonal Design. In general, the format of the orthogonal design is comparatively fixed if the number and the levels of the design variables have been given. That is, the optimal number of samples and its corresponding experimental factors can be determined if the ranges of the design and optimal levels are provided in advance. The distinctive features of schemes that are planned by orthogonal design are orthogonality and even distribution. Furthermore, orthogonal design requires fewer experiments in establishing a surrogate model, which improves the efficiency.(2) Central Composite Design (CCD). A complete scheme that is obtained by central composite design has three different experimental patterns. For example, if there are three design variables, the distribution of the sample points that are selected by CCD is illustrated in Figure 3. The eight points that are located at the summit of the cube are of the same experimental format, the six sample points that are located at the end of the center axis are of the same type, and the center point in red is a sample that needs to be repeated six times. Thus, for three variables, the optimal number of samples is twenty.Figure 3
Distribution of samples planned by CCD.


### 4.1.3. Establishment of Surrogate Model

Surrogate modeling has been broadly applied to many fields for its simple theory and excellent behavior. In the general case, the response surface methodology (RSM) is preferred by researchers for its simplicity and validity. Particularly, the RSM is more effective when dealing with multiobjective problems with multiple variables. Therefore, considering the characteristics of the proposed algorithm, the establishment of a linear, quadratic, or cubic polynomial surface is favored for describing the relationships between objective responses and design variables.Generally, a linear polynomial model is defined as follows:(25)y^=a+b1x1+b2x2+⋯+bmxm,where y^ is an approximate expression of objective or constraint functions, which is known as the response surface; xi is design variable i, with i=1,2,…,m; and m is the number of variables. The constant coefficients a,b1,b2,…,bm in (25) can be calculated in terms of the principle of least squares. Similarly, a quadratic polynomial model is usually constructed as in expression (26) and the cubic polynomial surface can be obtained in the same manner.(26)y^=a+∑1mbixi+∑1mbiixi2+∑i<kbikxixk.For simplicity, expression (26) can be converted to a linear form as follows:(27)y^=a+∑i=1WciXi,where W=m(m+3)/2 and X=(X1,X2,…,Xm) denote items with x1,x2,…,xm,x12,x22,…,xm2 and x1x2,x1x3,…,xm-1xm, respectively.To accurately estimate the precision of the approximated models, the correlation coefficientR2 and adjusted index Radj2 need to be calculated according to formulas (28) and (29), respectively.(28)R2=∑i=1My^i-y¯i2∑i=1Myi-y¯i2(29)Radj2=1-∑i=1Myi-y^i2M-1∑i=1Myi-y¯i2M-k-1,where M denotes the number of factors and y^i, y¯i  are the predictive value and the actual value, respectively, of the average measurement of y. The closer R2 and Radj2 are to 1, the better the surrogate model approximates the original problem. If they do not meet the accuracy requirement, the experimental points need to be reset to establish another model until the precision is satisfied.


### 4.1.4. Pareto Solutions Obtained by CBQNA

In this process, the approximate mathematical expressions of the primal multiobjective problems are first translated into a standard normalized multiobjective problem. Next, the Pareto-optimal solution set of the normalized problem is directly obtained by the CBQNA algorithm following the steps that are described in Section2. In practical projects, each Pareto solution corresponds to a design scheme, which can provide references for designers to optimize the original design.Based on these main steps, an effective optimization strategy for solving multiobjective optimization problems in engineering projects is in Figure4.Figure 4
Effective optimization strategy based on CBQNA.


## 4.1.1. Establishment of the Mathematical Model

To effectively solve optimization design problems, the optimization objectives and design variables need to be determined correctly. The preferred targets are often related to properties of the optimized structures. Analysis of sensitivity is usually applied when choosing decision variables, which greatly affects the targets. The design domains and constraint conditions are always determined according to the limiting dimension and strength requirements of the structure, scheme, or system. Hence, an accurate mathematical expression for the multiobjective optimization problem should be constructed based on the aforementioned analysis to solve the MOPs.


## 4.1.2. Design of Experiments

The establishment of an agent model with high precision is based on a good design scheme. The setting of reasonable factor levels and variable span play a significant role in the experimental design, especially in multifactor and multiobjective problems. Common methods for experimental design are orthogonal design and central composite design (CCD). In addition, Taguchi robust design has been applied in experimental design for its better performance [25]. To describe how to construct sample points and determine the optimal number of samples by orthogonal design and CCD, the procedural details of the two design methods can be separately summarized as follows.(1) Orthogonal Design. In general, the format of the orthogonal design is comparatively fixed if the number and the levels of the design variables have been given. That is, the optimal number of samples and its corresponding experimental factors can be determined if the ranges of the design and optimal levels are provided in advance. The distinctive features of schemes that are planned by orthogonal design are orthogonality and even distribution. Furthermore, orthogonal design requires fewer experiments in establishing a surrogate model, which improves the efficiency.(2) Central Composite Design (CCD). A complete scheme that is obtained by central composite design has three different experimental patterns. For example, if there are three design variables, the distribution of the sample points that are selected by CCD is illustrated in Figure 3. The eight points that are located at the summit of the cube are of the same experimental format, the six sample points that are located at the end of the center axis are of the same type, and the center point in red is a sample that needs to be repeated six times. Thus, for three variables, the optimal number of samples is twenty.Figure 3
Distribution of samples planned by CCD.


## 4.1.3. Establishment of Surrogate Model

Surrogate modeling has been broadly applied to many fields for its simple theory and excellent behavior. In the general case, the response surface methodology (RSM) is preferred by researchers for its simplicity and validity. Particularly, the RSM is more effective when dealing with multiobjective problems with multiple variables. Therefore, considering the characteristics of the proposed algorithm, the establishment of a linear, quadratic, or cubic polynomial surface is favored for describing the relationships between objective responses and design variables.Generally, a linear polynomial model is defined as follows:(25)y^=a+b1x1+b2x2+⋯+bmxm,where y^ is an approximate expression of objective or constraint functions, which is known as the response surface; xi is design variable i, with i=1,2,…,m; and m is the number of variables. The constant coefficients a,b1,b2,…,bm in (25) can be calculated in terms of the principle of least squares. Similarly, a quadratic polynomial model is usually constructed as in expression (26) and the cubic polynomial surface can be obtained in the same manner.(26)y^=a+∑1mbixi+∑1mbiixi2+∑i<kbikxixk.For simplicity, expression (26) can be converted to a linear form as follows:(27)y^=a+∑i=1WciXi,where W=m(m+3)/2 and X=(X1,X2,…,Xm) denote items with x1,x2,…,xm,x12,x22,…,xm2 and x1x2,x1x3,…,xm-1xm, respectively.To accurately estimate the precision of the approximated models, the correlation coefficientR2 and adjusted index Radj2 need to be calculated according to formulas (28) and (29), respectively.(28)R2=∑i=1My^i-y¯i2∑i=1Myi-y¯i2(29)Radj2=1-∑i=1Myi-y^i2M-1∑i=1Myi-y¯i2M-k-1,where M denotes the number of factors and y^i, y¯i  are the predictive value and the actual value, respectively, of the average measurement of y. The closer R2 and Radj2 are to 1, the better the surrogate model approximates the original problem. If they do not meet the accuracy requirement, the experimental points need to be reset to establish another model until the precision is satisfied.


## 4.1.4. Pareto Solutions Obtained by CBQNA

In this process, the approximate mathematical expressions of the primal multiobjective problems are first translated into a standard normalized multiobjective problem. Next, the Pareto-optimal solution set of the normalized problem is directly obtained by the CBQNA algorithm following the steps that are described in Section2. In practical projects, each Pareto solution corresponds to a design scheme, which can provide references for designers to optimize the original design.Based on these main steps, an effective optimization strategy for solving multiobjective optimization problems in engineering projects is in Figure4.Figure 4
Effective optimization strategy based on CBQNA.


## 4.2. Crashworthiness Optimization of Crash Box

With the boom of the automobile industry, the importance of passive safety has attracted increasing attention from researchers. To improve the energy absorption abilities of the main absorbers, such as the crash box, the longitudinal beam and door structure are vital to the protection of occupants. Thus, the analysis of the crashworthiness of the automobile energy absorbing box is essential and crucial.In this section, the proposed optimization strategy is applied to crashworthiness optimization of an automobile crash box. Consider a crash box that is composed of four thin-walled plates, as shown in Figure5. This thin-walled structure has been optimized by Yang et al. in study [32], but the impact condition that was used cannot correctly simulate the actual operation. Thus, a better simulation model under frontal collision is built in this paper.Figure 5
Bumper-crash box structure of the automobile.To effectively simulate the performance of the crash box under a low-speed impact, one end of each crash box has been fully constrained, while the other end is subject to front collision by a rigid wall of infinite mass at a low speed of 4 m/s along thex-axis, as shown in Figure 5. The original design of the crash is shown in Figure 6. The total impact time is 0.05 s.Figure 6
Original design of the crash box.


### 4.2.1. Mathematical Model for the MOP of the Crash Box Structure

The energy absorption and the maximal collision force, as the two major indices of the crashworthiness ability of the crash box, have been considered simultaneously. The crash box is a thin-walled structure, whose properties are affected dramatically by its wall thickness. In addition, when designing thin-walled structures, especially for automobile absorbers, the maximal crushing force and the mass should each be constrained below a maximum value. Thus, improving energy absorptionEin and reducing maximal impact force Fmax have been selected as simultaneous objectives, while the thicknesses of the four plates are chosen as design variables. Moreover, the mass of the entire structure and the maximal impact force are constrained by upper limits. Consequently, an MOP of the crash box that is based on crashworthiness is mathematically expressed as(30)min⁡F1=1Einx,F2=Fmaxxs.t.Fim≤F0;Mt≤M0⌻x∈xL,xU,where F1(x) and F2(x) are two objective functions of the vector x=(t1,t2,t3,t4), which consists of four decision variables, which are subject to boundary condition ti∈[1,3]mm,i=1,2,3,4. The conditions Fmax≤F0 and Mt≤M0 constrain the impact force and the total mass should be maintained below an upper value. In this paper, we set F0=174 kN and M0=4.2 kg.


### 4.2.2. Experiment Design

The thicknesses of the four thin-walled plates are concurrently considered by DOE. The combined orthogonal experimental method is utilized to construct a 17-group experimental dataset in the design domain.


### 4.2.3. Construction of Response Surface Models

The response surface methodology is utilized to construct surrogate models instead of the original complicated problem. The mathematical expressions for the two objectives in terms of the four variables are stated as (31) and (32). In addition, the approximate mode of the constraint function is presented in expression (33).(31)Einx=-3397.09+2967.8t1+3368.16t2+1300.84t3+254.58t4-157.9t12-211.4t22-31.41t32-38.69t42-684.18t1t2-93.37t1t3-307.1t2t3+85.12t2t4-65.4t3t4-38.7t1t4;(32)Fmaxx=-76.71-13.62t1-10.5t2+107.4t3+142.8t4+28.14t12+24.44t22-18.43t32-18.38t42-11.13t1t2-25.40t1t3-6.38t2t3-39.58t2t4+13.99t3t4-10.96t1t4;(33)Mtx=-0.0128+0.389t1+1.19t2+0.346t3+0.183t4.After calculation, the coefficients of determination of energy absorption, maximal impact force, and mass constraint are 98.60%, 99.11%, and 1, respectively, while the adjusted coefficients of determination are 95.6%, 97.11%, and 1, respectively. The calculated indicators demonstrate that the approximate models of the two objective and constraint functions have excellent precision.


### 4.2.4. MOP by CBQNA

The hybrid algorithm CBQNA, which is introduced in Section2, is applied to solve the optimization design problem. For adapting the surrogate model to the novel algorithm, a normalization process is implemented first. Thus, problem (30) can be converted into a standard mathematical problem, which is expressed as follows:(34)minF1x=EmaxEinxF2x=FmaxxFas.t.1mm≤ti≤3mm,i=1,2,3,4.Fmaxx≤174kN;Mtx≤4.2kg,where Emax and Fa are the maximal value of the energy absorption and average value of the maximal collision force, respectively.After optimization, 100 groups of solutions of the optimization problem are obtained, which takes only 5.6 seconds on average. The Pareto-optimal front, which is composed of 100 Pareto points, is presented in Figure7, and the relationship between energy absorption and maximal impact force is depicted in Figure 8.Figure 7
Pareto-optimal front obtained by CBQNA.Figure 8
Relationship between the energy absorption and the maximal impact force.In terms of the properties of energy absorption ability and maximal collision force, the optimal design is specified byt1=2.0,t2=1.98,t3=1.98, and t4=1.98. Regarding the lightweight design, the mass of the thin-walled plates is also a significant factor in selecting the best design. The relationship between the energy absorption, maximal impact force, and total mass of the bumper absorbing box is depicted in Figure 9. The comparison between the original design and the optimal result in terms of energy absorption, maximal collision force, and mass is presented in Table 2.Table 2
Comparison of design variables and the performances before and after optimization.
t
i
/
m
m
E
i
n
/
J
F
m
a
x
/
k
N
M
a
/
k
g
t
1
t
2
t
3
t
4
Before
1.4
1.8
2.3
1.8
5605.2
173.25
3.83
After
1.86
1.78
2.0
1.96
5777.8
146.64
3.89Figure 9
Objective responses and total mass of the crash box structure.Table2 shows that, compared with the original design, the performance improves greatly after optimization. The energy absorption improves by 3.08% and the maximal impact force reduces by 14.87%, while the total mass increases slightly by 1.57%. The collision conditions of the crash box at 0.01 s, 0.02 s, 0.03 s, and 0.04 s are presented in Figure 10. From the collision deformation, the optimized crash box structure achieves an ideal deformation and satisfies the practical engineering requirements.Figure 10
Deformation of the crash box at different times.


## 4.2.1. Mathematical Model for the MOP of the Crash Box Structure

The energy absorption and the maximal collision force, as the two major indices of the crashworthiness ability of the crash box, have been considered simultaneously. The crash box is a thin-walled structure, whose properties are affected dramatically by its wall thickness. In addition, when designing thin-walled structures, especially for automobile absorbers, the maximal crushing force and the mass should each be constrained below a maximum value. Thus, improving energy absorptionEin and reducing maximal impact force Fmax have been selected as simultaneous objectives, while the thicknesses of the four plates are chosen as design variables. Moreover, the mass of the entire structure and the maximal impact force are constrained by upper limits. Consequently, an MOP of the crash box that is based on crashworthiness is mathematically expressed as(30)min⁡F1=1Einx,F2=Fmaxxs.t.Fim≤F0;Mt≤M0⌻x∈xL,xU,where F1(x) and F2(x) are two objective functions of the vector x=(t1,t2,t3,t4), which consists of four decision variables, which are subject to boundary condition ti∈[1,3]mm,i=1,2,3,4. The conditions Fmax≤F0 and Mt≤M0 constrain the impact force and the total mass should be maintained below an upper value. In this paper, we set F0=174 kN and M0=4.2 kg.


## 4.2.2. Experiment Design

The thicknesses of the four thin-walled plates are concurrently considered by DOE. The combined orthogonal experimental method is utilized to construct a 17-group experimental dataset in the design domain.


## 4.2.3. Construction of Response Surface Models

The response surface methodology is utilized to construct surrogate models instead of the original complicated problem. The mathematical expressions for the two objectives in terms of the four variables are stated as (31) and (32). In addition, the approximate mode of the constraint function is presented in expression (33).(31)Einx=-3397.09+2967.8t1+3368.16t2+1300.84t3+254.58t4-157.9t12-211.4t22-31.41t32-38.69t42-684.18t1t2-93.37t1t3-307.1t2t3+85.12t2t4-65.4t3t4-38.7t1t4;(32)Fmaxx=-76.71-13.62t1-10.5t2+107.4t3+142.8t4+28.14t12+24.44t22-18.43t32-18.38t42-11.13t1t2-25.40t1t3-6.38t2t3-39.58t2t4+13.99t3t4-10.96t1t4;(33)Mtx=-0.0128+0.389t1+1.19t2+0.346t3+0.183t4.After calculation, the coefficients of determination of energy absorption, maximal impact force, and mass constraint are 98.60%, 99.11%, and 1, respectively, while the adjusted coefficients of determination are 95.6%, 97.11%, and 1, respectively. The calculated indicators demonstrate that the approximate models of the two objective and constraint functions have excellent precision.


## 4.2.4. MOP by CBQNA

The hybrid algorithm CBQNA, which is introduced in Section2, is applied to solve the optimization design problem. For adapting the surrogate model to the novel algorithm, a normalization process is implemented first. Thus, problem (30) can be converted into a standard mathematical problem, which is expressed as follows:(34)minF1x=EmaxEinxF2x=FmaxxFas.t.1mm≤ti≤3mm,i=1,2,3,4.Fmaxx≤174kN;Mtx≤4.2kg,where Emax and Fa are the maximal value of the energy absorption and average value of the maximal collision force, respectively.After optimization, 100 groups of solutions of the optimization problem are obtained, which takes only 5.6 seconds on average. The Pareto-optimal front, which is composed of 100 Pareto points, is presented in Figure7, and the relationship between energy absorption and maximal impact force is depicted in Figure 8.Figure 7
Pareto-optimal front obtained by CBQNA.Figure 8
Relationship between the energy absorption and the maximal impact force.In terms of the properties of energy absorption ability and maximal collision force, the optimal design is specified byt1=2.0,t2=1.98,t3=1.98, and t4=1.98. Regarding the lightweight design, the mass of the thin-walled plates is also a significant factor in selecting the best design. The relationship between the energy absorption, maximal impact force, and total mass of the bumper absorbing box is depicted in Figure 9. The comparison between the original design and the optimal result in terms of energy absorption, maximal collision force, and mass is presented in Table 2.Table 2
Comparison of design variables and the performances before and after optimization.
t
i
/
m
m
E
i
n
/
J
F
m
a
x
/
k
N
M
a
/
k
g
t
1
t
2
t
3
t
4
Before
1.4
1.8
2.3
1.8
5605.2
173.25
3.83
After
1.86
1.78
2.0
1.96
5777.8
146.64
3.89Figure 9
Objective responses and total mass of the crash box structure.Table2 shows that, compared with the original design, the performance improves greatly after optimization. The energy absorption improves by 3.08% and the maximal impact force reduces by 14.87%, while the total mass increases slightly by 1.57%. The collision conditions of the crash box at 0.01 s, 0.02 s, 0.03 s, and 0.04 s are presented in Figure 10. From the collision deformation, the optimized crash box structure achieves an ideal deformation and satisfies the practical engineering requirements.Figure 10
Deformation of the crash box at different times.


## 4.3. Optimization Design of an S-Shaped Thin-Walled Beam

The frontal rail that is connected with the crash box is also a main deformation and energy absorption structure in a frontal crash, which is shown in Figure11.Figure 11
Simplified model of an S-shaped front rail structure.Generally, the vehicle front rail is a thin-walled structure whose deformation and energy absorption characteristics are similar to those of a general rectangular thin-walled beam. Thus, the crashworthiness of the thin-walled beam has dramatic significance to the frontal crash safety of the automobile. In this paper, an attempt has been made to study an S-shaped simplified model under frontal impact conditions.


### 4.3.1. Establishment of the Finite-Element Model

Firstly, the geometry model of S-shaped structure is established with Catia, then its finite-element model is simulated by Hypermesh, and the simulation process is carried by LS-DYNA. The FE model of the S-shaped thin-walled beam with rectangular cross-section is displayed in Figure12.The S-shaped beam is impacted by a rigid wall that weights 450 kg at the speed of 10 m/s. As shown in Figure 13, the rear end of the beam is fully constrained and the other end is allowed to move only along the X-direction. The material that is utilized for the S-shaped beam is mild steel with the following material properties: Young’s modulus = 206 Gpa, Poisson’s ratio = 0.3, and density = 7800 kg/m3. The relationship between stress and strain rate (shown in Figure 14) is also considered.Figure 12
FE model of the S-shaped beam.Figure 13
Boundary condition for S-shaped front rail.Figure 14
Typical engineering stress-strain curve of the material [47].


### 4.3.2. Mathematical Model for MOP of the Front Rail

The cross-sectional dimensions of the thin-walled structure are the main factor that affects the properties of energy absorption and maximal impact resistance. The widthw, height h, and thickness t of the cross-section have been chosen as design variables in this study, and the cross-section of an S-shaped thin-walled structure is shown in Figure 15.Figure 15
Cross-section of an S-shaped beam.In addition, energy absorption and maximal impact force are major evaluation indices of the crashworthiness ability of the S-shaped front rail beam.Hence, the mathematical model for multiobjective optimization of the front rail, which consists of two objective functions in terms of three design variables, is constructed as(35)min⁡Fx=F1x=1ExF2x=Fxs.t.hmin≤h≤hmaxwmin≤w≤wmaxtmin≤t≤tmax,where x=(w,h,t)T is the vector of design variables; F(x) is a vector function that consists of two objectives, namely, F1(x) and F2(x); hmin and hmax denote the lower and upper limits of the height of the cross-section; wmin and wmax denote the lower and upper limits of the width; and tmin and tmax are the lower and upper limits of the thickness of the wall, respectively. In this study, to satisfy the limits of the structure dimensions and the strength demands, hmin=136mm,hmax=146mm,wmin=56mm,wmax=68mm,tmin=1.4mm, and tmax=2.2mm.


### 4.3.3. Design of Experiments

Following the novel optimization strategy’s procedures, orthogonal design is utilized to generate 25 design sampling points (i.e., four levels for each of the three design variables), which are clearly presented in Table3 with their corresponding responses.Table 3
Data for 25 groups, obtained by orthogonal experimental design.
Number
Factors
E  (J)
F  (kN)
w
h
t
(1)
58
136
1.4
11444.7
182.067
(2)
58
139
1.6
12943.7
214.618
(3)
58
142
1.8
14335.4
235.309
(4)
58
145
2.0
15644
251.437
(5)
58
148
2.2
16845.9
287.804
(6)
60
136
1.6
12910.2
212.551
(7)
60
139
1.8
14301
233.593
(8)
60
142
2.0
15603.6
248.585
(9)
60
145
2.2
16805.2
272.659
(10)
60
148
1.4
12257.2
196.856
(11)
62
136
1.8
14242.8
221.679
(12)
62
139
2.0
15547.9
246.477
(13)
62
142
2.2
16643.2
299.146
(14)
62
145
1.4
12236.6
197.617
(15)
62
148
1.6
13765.4
231.555
(16)
64
136
2.0
15383.9
259.826
(17)
64
139
2.2
16578.2
299.329
(18)
64
142
1.4
12200.4
194.745
(19)
64
145
1.6
13730.5
229.025
(20)
64
148
1.8
15148.8
241.987
(21)
66
136
2.2
16513.7
299.543
(22)
66
139
1.4
12157.1
195.991
(23)
66
142
1.6
13679.3
224.678
(24)
66
145
1.8
15094.3
234.490
(25)
66
148
2.0
16400.2
272.619


### 4.3.4. Agent Model and Precision Analysis

A surrogate model of two objectives for crashworthiness optimization for the frontal rail is also established to lower the complexity of the original problem. According to the 25 plans of the DOE, which are shown in Table3, the quadratic response surface models of the energy absorption and maximal impact force responses were constructed in expressions (36) and (37), respectively.(36)E′x=-70973.4+958.656w+614.173h-3.827wh+57.023wt+63.274ht-3.582w2-1.4842h2-1890.174t2(37)F′x=111.854-20.210w+7.560h-0.025wh+2.731wt-1.479ht+0.165w2-0.009h2+43.110t2.The surrogate model’s fitness was evaluated based on coefficientR2 and its adjusted form Radj2. The values of coefficient R2 for the energy absorption and collision force reach 99% and 98%, respectively, and the correction coefficient Radj2 values are 99% and 95%, respectively. The result indicates that the approximation ability of this surrogate model is very good.


### 4.3.5. Optimization Design

To improve the efficiency of optimization, the hybrid algorithm CBQNA is utilized to minimize the two response functions simultaneously. First, based on the new algorithm, a normalization procedure is applied to unify the magnitudes of the two responses functions. Thus, problem (35) is converted into a standard biobjective optimization problem with inequality constraints, which is expressed as follows:(38)min⁡Fx=F1x=EaE′xF2x=F′xFa(39)s.t.h-136146-h≥0;w-5866-w≥0,t-1.42.2-t≥0;where Ea=14496.5J and Fa=234.9 kN are the average value of energy absorption and the maximal collision force for the 25 groups, respectively.


### 4.3.6. Optimization Results Solved by CBQNA

Problem (38) is solved by the new proposed algorithm with random initial value in accordance with constraints. To maintain the diversity of design selection, 200 groups of optimal solutions are acquired regularly by adjusting weighting factors. The Pareto-optimal front, which is constructed with 200 solutions, is depicted in Figure 16. In addition, the Pareto-optimal front of the front rail collision problem is shown in Figure 17. When determining the dimensions of the front rail, we should focus on not only the crashworthiness ability but also the weight requirement of the automobile. Thus, the mass of the thin-walled beam is also considered in evaluating the crashworthiness performance in Table 4.Table 4
Optimal solutions and their corresponding objective responses.
Design variables
E
(
J
)
F
(
k
N
)
M
(
k
g
)
w
h
t
(1)
63.18507
145.4626
2.1208
17081.0
272.889
7.19
(2)
63.23390
145.5450
1.8717
15538.6
245.962
6.36
(3)
63.23933
145.5548
2.0046
16087.7
263.545
6.80
(4)
63.23917
145.5540
1.9672
16087.7
263.545
6.67
(5)
63.17605
145.4450
2.1677
17048.5
267.451
7.36
(6)
63.05213
145.2380
1.8668
15485.6
245.663
6.32
(7)
63.10070
145.3174
1.7284
14311.9
233.394
5.89
(8)
63.16113
145.4234
1.8880
15506.7
245.863
6.43
(9)
63.22119
145.5222
1.9085
15511.4
245.155
6.46
(10)
63.21939
145.5207
1.8563
15509.1
243.437
6.33
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
(71)
63.20308
145.4854
1.8738
15506.2
245.319
6.38
(72)
63.14951
145.3912
1.7064
14323.5
233.614
5.81
(73)
63.18035
145.4439
2.1700
17036.8
276.266
7.38
(74)
63.07671
145.2661
1.9681
16030.8
257.892
6.68
(75)
63.18094
145.4458
1.9650
16049.1
258.247
6.69
(76)
63.14457
145.3848
2.0132
16037.5
260.135
7.02
(77)
63.19502
145.4693
2.1381
16586.5
266.359
7.26
(78)
63.20000
145.4776
1.8571
15677.5
240.558
6.23
(79)
63.11433
145.3301
2.1676
17037.7
265.008
7.35
(80)
63.20090
145.4790
1.9712
16056.5
257.411
6.70
(81)
63.19128
145.4618
2.0219
16054.8
264.999
6.86
(82)
63.15801
145.4045
1.9299
15504.9
246.552
6.56
(83)
63.12779
145.3543
1.9100
15507.3
245.189
6.49
(84)
63.20623
145.4870
1.9823
16052.2
258.096
6.73
(85)
63.18679
145.4535
1.9746
16055.6
257.415
6.70
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…Figure 16
Pareto-optimal front obtained by CBQNA.Figure 17
Relationship between the energy absorption and the maximal impact force.According to the properties of the best plans, the cross-sectional dimensions are selected, and the thickness of the wall is set tox=63.2,146,1.8mm, as shown in line 78 in Table 4. After the simulation by software LS-DYNA, the energy absorption and collision force responses, the specific energy absorption (SEA), and the mass of the structure, compared with the original design, are presented in Table 5. Figure 18 illustrates the deformation of the S-shaped frontal rail thin-walled structure changes over time.Table 5
Design variables and response values before and after optimization.
Design variables (mm)
E
(
J
)
F
(
k
N
)
M
(
k
g
)
w
h
t
Original Design
62.0
142
2.00
15616.8
247.9
6.45
Optimized Result
63.2
146
1.80
15077.5
240.5
6.03Figure 18
Deformation of the S-shaped beam changes over time.Although, compared with the original design, the internal energy absorption of the optimal design is decreased slightly, the maximal impact and the mass are decreased by 3.0% and 6.5%, respectively. Moreover, the SEA, which is the main indicator of the energy absorption ability, drops by 3%, which signifies that the crashworthiness of the frontal rail is improved. Therefore, all performance indicators of the optimal result meet the requirements.


## 4.3.1. Establishment of the Finite-Element Model

Firstly, the geometry model of S-shaped structure is established with Catia, then its finite-element model is simulated by Hypermesh, and the simulation process is carried by LS-DYNA. The FE model of the S-shaped thin-walled beam with rectangular cross-section is displayed in Figure12.The S-shaped beam is impacted by a rigid wall that weights 450 kg at the speed of 10 m/s. As shown in Figure 13, the rear end of the beam is fully constrained and the other end is allowed to move only along the X-direction. The material that is utilized for the S-shaped beam is mild steel with the following material properties: Young’s modulus = 206 Gpa, Poisson’s ratio = 0.3, and density = 7800 kg/m3. The relationship between stress and strain rate (shown in Figure 14) is also considered.Figure 12
FE model of the S-shaped beam.Figure 13
Boundary condition for S-shaped front rail.Figure 14
Typical engineering stress-strain curve of the material [47].


## 4.3.2. Mathematical Model for MOP of the Front Rail

The cross-sectional dimensions of the thin-walled structure are the main factor that affects the properties of energy absorption and maximal impact resistance. The widthw, height h, and thickness t of the cross-section have been chosen as design variables in this study, and the cross-section of an S-shaped thin-walled structure is shown in Figure 15.Figure 15
Cross-section of an S-shaped beam.In addition, energy absorption and maximal impact force are major evaluation indices of the crashworthiness ability of the S-shaped front rail beam.Hence, the mathematical model for multiobjective optimization of the front rail, which consists of two objective functions in terms of three design variables, is constructed as(35)min⁡Fx=F1x=1ExF2x=Fxs.t.hmin≤h≤hmaxwmin≤w≤wmaxtmin≤t≤tmax,where x=(w,h,t)T is the vector of design variables; F(x) is a vector function that consists of two objectives, namely, F1(x) and F2(x); hmin and hmax denote the lower and upper limits of the height of the cross-section; wmin and wmax denote the lower and upper limits of the width; and tmin and tmax are the lower and upper limits of the thickness of the wall, respectively. In this study, to satisfy the limits of the structure dimensions and the strength demands, hmin=136mm,hmax=146mm,wmin=56mm,wmax=68mm,tmin=1.4mm, and tmax=2.2mm.


## 4.3.3. Design of Experiments

Following the novel optimization strategy’s procedures, orthogonal design is utilized to generate 25 design sampling points (i.e., four levels for each of the three design variables), which are clearly presented in Table3 with their corresponding responses.Table 3
Data for 25 groups, obtained by orthogonal experimental design.
Number
Factors
E  (J)
F  (kN)
w
h
t
(1)
58
136
1.4
11444.7
182.067
(2)
58
139
1.6
12943.7
214.618
(3)
58
142
1.8
14335.4
235.309
(4)
58
145
2.0
15644
251.437
(5)
58
148
2.2
16845.9
287.804
(6)
60
136
1.6
12910.2
212.551
(7)
60
139
1.8
14301
233.593
(8)
60
142
2.0
15603.6
248.585
(9)
60
145
2.2
16805.2
272.659
(10)
60
148
1.4
12257.2
196.856
(11)
62
136
1.8
14242.8
221.679
(12)
62
139
2.0
15547.9
246.477
(13)
62
142
2.2
16643.2
299.146
(14)
62
145
1.4
12236.6
197.617
(15)
62
148
1.6
13765.4
231.555
(16)
64
136
2.0
15383.9
259.826
(17)
64
139
2.2
16578.2
299.329
(18)
64
142
1.4
12200.4
194.745
(19)
64
145
1.6
13730.5
229.025
(20)
64
148
1.8
15148.8
241.987
(21)
66
136
2.2
16513.7
299.543
(22)
66
139
1.4
12157.1
195.991
(23)
66
142
1.6
13679.3
224.678
(24)
66
145
1.8
15094.3
234.490
(25)
66
148
2.0
16400.2
272.619


## 4.3.4. Agent Model and Precision Analysis

A surrogate model of two objectives for crashworthiness optimization for the frontal rail is also established to lower the complexity of the original problem. According to the 25 plans of the DOE, which are shown in Table3, the quadratic response surface models of the energy absorption and maximal impact force responses were constructed in expressions (36) and (37), respectively.(36)E′x=-70973.4+958.656w+614.173h-3.827wh+57.023wt+63.274ht-3.582w2-1.4842h2-1890.174t2(37)F′x=111.854-20.210w+7.560h-0.025wh+2.731wt-1.479ht+0.165w2-0.009h2+43.110t2.The surrogate model’s fitness was evaluated based on coefficientR2 and its adjusted form Radj2. The values of coefficient R2 for the energy absorption and collision force reach 99% and 98%, respectively, and the correction coefficient Radj2 values are 99% and 95%, respectively. The result indicates that the approximation ability of this surrogate model is very good.


## 4.3.5. Optimization Design

To improve the efficiency of optimization, the hybrid algorithm CBQNA is utilized to minimize the two response functions simultaneously. First, based on the new algorithm, a normalization procedure is applied to unify the magnitudes of the two responses functions. Thus, problem (35) is converted into a standard biobjective optimization problem with inequality constraints, which is expressed as follows:(38)min⁡Fx=F1x=EaE′xF2x=F′xFa(39)s.t.h-136146-h≥0;w-5866-w≥0,t-1.42.2-t≥0;where Ea=14496.5J and Fa=234.9 kN are the average value of energy absorption and the maximal collision force for the 25 groups, respectively.


## 4.3.6. Optimization Results Solved by CBQNA

Problem (38) is solved by the new proposed algorithm with random initial value in accordance with constraints. To maintain the diversity of design selection, 200 groups of optimal solutions are acquired regularly by adjusting weighting factors. The Pareto-optimal front, which is constructed with 200 solutions, is depicted in Figure 16. In addition, the Pareto-optimal front of the front rail collision problem is shown in Figure 17. When determining the dimensions of the front rail, we should focus on not only the crashworthiness ability but also the weight requirement of the automobile. Thus, the mass of the thin-walled beam is also considered in evaluating the crashworthiness performance in Table 4.Table 4
Optimal solutions and their corresponding objective responses.
Design variables
E
(
J
)
F
(
k
N
)
M
(
k
g
)
w
h
t
(1)
63.18507
145.4626
2.1208
17081.0
272.889
7.19
(2)
63.23390
145.5450
1.8717
15538.6
245.962
6.36
(3)
63.23933
145.5548
2.0046
16087.7
263.545
6.80
(4)
63.23917
145.5540
1.9672
16087.7
263.545
6.67
(5)
63.17605
145.4450
2.1677
17048.5
267.451
7.36
(6)
63.05213
145.2380
1.8668
15485.6
245.663
6.32
(7)
63.10070
145.3174
1.7284
14311.9
233.394
5.89
(8)
63.16113
145.4234
1.8880
15506.7
245.863
6.43
(9)
63.22119
145.5222
1.9085
15511.4
245.155
6.46
(10)
63.21939
145.5207
1.8563
15509.1
243.437
6.33
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
(71)
63.20308
145.4854
1.8738
15506.2
245.319
6.38
(72)
63.14951
145.3912
1.7064
14323.5
233.614
5.81
(73)
63.18035
145.4439
2.1700
17036.8
276.266
7.38
(74)
63.07671
145.2661
1.9681
16030.8
257.892
6.68
(75)
63.18094
145.4458
1.9650
16049.1
258.247
6.69
(76)
63.14457
145.3848
2.0132
16037.5
260.135
7.02
(77)
63.19502
145.4693
2.1381
16586.5
266.359
7.26
(78)
63.20000
145.4776
1.8571
15677.5
240.558
6.23
(79)
63.11433
145.3301
2.1676
17037.7
265.008
7.35
(80)
63.20090
145.4790
1.9712
16056.5
257.411
6.70
(81)
63.19128
145.4618
2.0219
16054.8
264.999
6.86
(82)
63.15801
145.4045
1.9299
15504.9
246.552
6.56
(83)
63.12779
145.3543
1.9100
15507.3
245.189
6.49
(84)
63.20623
145.4870
1.9823
16052.2
258.096
6.73
(85)
63.18679
145.4535
1.9746
16055.6
257.415
6.70
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…
…Figure 16
Pareto-optimal front obtained by CBQNA.Figure 17
Relationship between the energy absorption and the maximal impact force.According to the properties of the best plans, the cross-sectional dimensions are selected, and the thickness of the wall is set tox=63.2,146,1.8mm, as shown in line 78 in Table 4. After the simulation by software LS-DYNA, the energy absorption and collision force responses, the specific energy absorption (SEA), and the mass of the structure, compared with the original design, are presented in Table 5. Figure 18 illustrates the deformation of the S-shaped frontal rail thin-walled structure changes over time.Table 5
Design variables and response values before and after optimization.
Design variables (mm)
E
(
J
)
F
(
k
N
)
M
(
k
g
)
w
h
t
Original Design
62.0
142
2.00
15616.8
247.9
6.45
Optimized Result
63.2
146
1.80
15077.5
240.5
6.03Figure 18
Deformation of the S-shaped beam changes over time.Although, compared with the original design, the internal energy absorption of the optimal design is decreased slightly, the maximal impact and the mass are decreased by 3.0% and 6.5%, respectively. Moreover, the SEA, which is the main indicator of the energy absorption ability, drops by 3%, which signifies that the crashworthiness of the frontal rail is improved. Therefore, all performance indicators of the optimal result meet the requirements.


## 4.4. Analysis of Optimal Results

The optimization strategy that is based on CBQNA has been applied to the optimization design of the crash box and the S-shaped front rail. The optimized results demonstrate that the optimization strategy that is based on the CBQNA algorithm is valuable in solving MOPs in engineering projects. To evaluate the efficiency of the proposed method, the two engineering examples in this paper are each executed 30 times. The calculation results show that it takes approximately 5.6 s to optimize the crash box by this novel algorithm, compared to approximately 198 s by NSWFA [32]. Furthermore, the constraint condition of the crash box that is simulated in this paper is more consistent with the actual frontal collision than that in [32]. The calculation time of the optimization of the S-shaped front rail is approximately 32.8 s for 200 groups’ optimal solutions, which demonstrates that the optimization strategy that is based on the proposed algorithm is effective and promising in solving complicated engineering problems in which the objective function’s second derivative is difficult to calculate or does not exist.


## 5. Conclusions

This article proposed a new hybrid algorithm, known as CBQNA, for solving multiobjective problems and a surrogate-modeling-based optimization strategy to address MOPs in engineering projects. The hybrid algorithm combines a weighted sum method and a cautious BFGS method to determine a descent search direction. Avoiding the calculation of the Hessian matrix of each objective function in every iteration dramatically reduces the time consumption. The algorithm CBQNA has excellent performance in terms of obtaining uniformly distributed solutions and high computational efficiency, which has been proven by numerical experiments. The optimization strategy has been applied to crashworthiness optimal design of an automobile crash box and an S-shaped front rail structure. The outstanding performance of the surrogate-modeling-based strategy that was proposed in this paper can be summarized as follows. First, because of the high calculation speed of the hybrid algorithm, the engineering problem can be solved within the allowed time, which is very significant in the actual design process. Second, the optimal result is superior to the original design in terms of crashworthiness ability, which is crucial to vehicle passive safety. In our future work, evolutionary algorithms will be combined with the proposed algorithm to generate better initial value for definitely avoiding sinking into pseudo-optimal solution and the convergence of the proposed hybrid algorithm will be theoretically proven.



---

*Source: 5316379-2018-01-22.xml*