# An Efficient Imperialist Competitive Algorithm for Solving the QFD Decision Problem

Authors: Xue Ji; Qi Gao; Fupeng Yin; Hengdong Guo
Journal: Mathematical Problems in Engineering
 (2016)

---

## Abstract

It is an important QFD decision problem to determine the engineering characteristics and their corresponding actual fulfillment levels. With the increasing complexity of actual engineering problems, the corresponding QFD matrixes become much huger, and the time spent on analyzing these matrixes and making decisions will be unacceptable. In this paper, a solution for efficiently solving the QFD decision problem is proposed. The QFD decision problem is reformulated as a mixed integer nonlinear programming (MINLP) model, which aims to maximize overall customer satisfaction with the consideration of the enterprises’ capability, cost, and resource constraints. And then an improved algorithm G-ICA, a combination of Imperialist Competitive Algorithm (ICA) and genetic algorithm (GA), is proposed to tackle this model. The G-ICA is compared with other mature algorithms by solving 7 numerical MINLP problems and 4 adapted QFD decision problems with different scales. The results verify a satisfied global optimization performance and time performance of the G-ICA. Meanwhile, the proposed algorithm’s better capabilities to guarantee decision-making accuracy and efficiency are also proved.


---

## Body

## 1. Introduction

Quality Function Deployment (QFD), as a well-known planning methodology originated in the late 1960s in Japan by Akao and Mazur [1], focuses on translating customer requirements (CRs) into a series of relevant design requirements. During this planning process, a lot of decisions should be made. One of them is to determine the engineering characteristics (ECs) to be implemented among numerous possible selections and their corresponding actual fulfillment levels, which not only directly influences the deployment of other design requirements and the final quality of the product but also influences overall customer satisfaction and the distribution of resources.In practice, it is a complex process to analyze and organize a lot of information related to the QFD decision problem, including multifacets of customer needs, market driven factors, price/revenue streams of enterprises, and various regulation compliance. Multiple variables and complex relationships bring many difficulties for decision makers to make decisions. Mathematical programming, an optimization technique, such as integer programming (IP), linear programming (LP), nonlinear programming (NLP), and mixed integer nonlinear programming (MINLP), has been widely used in helping people solve these difficulties. Firstly, in order to clarify those complex relationships, the problem is abstracted as a mathematical model. Then the model is calculated by an algorithm to offer a result for decision makers as a valuable reference [2–13].However, a rapidly changing decision environment makes people not only satisfied with obtaining a decision outcome but also want to seek for a higher decision-making efficiency, in other words, a shorter time spent on making decisions. Especially in such a society with highly developed economy, technology, and culture, customer requirements tend to be more and more diverse, product systems become gradually complex, and science technologies are continuously developing. All of these may lead to a huge and complex QFD matrix. As an important step in tackling the above QFD decision problem, the selection of algorithms influences not only the accuracy of decisions but also total decision time. And if the calculating time of the algorithm is too long, the decision efficiency and time cannot be guaranteed.The QFD decision problem in this paper is reformulated as a MINLP model aiming to maximize overall customer satisfaction with the consideration of the enterprises’ capability, cost, and resource constraints. Actually, many decision problems in various application fields, such as chemical engineering [14, 15], gas/water distribution networks optimal design [16, 17], and engineering design [18], can be formulated as MINLP models. MINLP is a NP-hard problem and the calculating algorithms to MINLPs have been developed since 1980s, which can be broadly divided into two main classes: deterministic algorithms and stochastic algorithms. Deterministic algorithms, such as branch and bound method, outer approximation method, and extended cutting plane method, can provide accurate solutions for problems in a continuous space but require an additional effort of analysis (such as the gradient analysis of objective function and constraints), which may considerably add to the work [19, 20]. Considering the existence of integer variables in QFD decision problems and easy-operation of engineers, stochastic algorithms, which are more capable of handling NP-hard problems and need not much professional mathematical knowledge, will be focused on in this paper. A lot of stochastic algorithms have been proposed. An algorithm called M-SIMPSA suitable for the optimization of MINLP is presented by Cardoso et al. [19]. Mohan and Nguyen [21] developed an improved computational algorithm called RST2ANU, which was based on the original RST2U [22]. Yiqing et al. [23] improved the particle swarm optimization (PSO) algorithm for solving nonconvex NLP/MINLP problem with equality and/or inequality constraints. Schlüter et al. [24] introduced two novel extensions for ant colony optimization (ACO) framework. Deep et al. [25] proposed a real coded genetic algorithm named MI-LXPM.In this paper, Imperialist Competitive Algorithm (ICA), one of the stochastic algorithms proposed by Atashpaz-Gargari and Lucas [26], which simulates the social-political process of imperialism and imperialist competition, is applied to solve MINLP models. The key features of ICA are its fast convergent rate to reach global optimum, which has been proved in dealing with various optimization problems, such as structure design [27, 28], distribution network optimization [29, 30], and process planning and scheduling [31]. The advantages of ICA are beneficial to improvement of decision efficiency. In this QFD decision problem, considering the existence of integer variables in this problem, an improvement by combining genetic algorithm (GA) is carried out because the original ICA is only applied to real variables. The improved ICA, called G-IGA in the following, is compared with other mature algorithms to MINLPs, including RST2U and RST2ANU mentioned above, by calculating 7 simple numerical MINLP problems and 4 different scale QFD decision problems, and the results verify the advantages of G-ICA.The remainder of the paper is organized as follows: the QFD decision problem is defined in Section2. The G-ICA to resolve the problem is introduced in Section 3. Some numerical experiments and case studies are shown in Section 4, and Section 5 summarizes the paper.


## 2. Problem Definition

The QFD matrix is just as in Figure1. Before the decision problem based on QFD is defined, some assumptions should be made in advance:(i)
The QFD matrix for a certain product, includingm CRs, n possible ECs, weight vectors, relationship, and correlationship matrixes, is already established.(ii)
Data about capability (lower and upper bounds of enterprise fulfillment levels for every EC), cost (maximum total budget and upper bound of investment cost for every EC), and resource (maximum total resources and required resource to implement every EC) is offered by the enterprise.(iii)
The fulfillment level is a percentage from 0 to 100%. For every EC, actual fulfillment level is only influenced by the corresponding investment cost, whose value is from 0 to 10 with no units. And the enterprise can determine the influence relationship based on long-term practical experience all by itself.(iv)
There are no natural units or maximum values for customer satisfaction. It is the result of synthesizing all the actual fulfillment levels of engineering characteristics to be implemented.Figure 1
The QFD matrix.Then the mathematical model is expressed as follows:(1)maxCS=∑j=1nNj·ρj·gj(2)subjecttogLj≤gj≤gUj(3)0≤cj≤cUj(4)∑j=1nNj·cj≤B(5)∑j=1nNj·pj≤P(6)Nj∈0,1(7)gj=gjcjin which the decision variables and parameters are defined as follows:CS: the overall customer satisfactionj: index for possible engineering characteristics, j=1,2,…,nNj: if engineering characteristic j is selected, Nj is taken as 1; otherwise, Nj is taken as 0ρj: the final importance weight of the engineering characteristic jgj: the actual fulfillment level for engineering characteristic j, gj∈[0,1]gLj: the lower bound of enterprise fulfillment level of engineering characteristic j, gLj∈[0,1]gUj: the upper bound of enterprise fulfillment level of engineering characteristic j, gUj∈[0,1]cj: the investment cost for engineering characteristic j, cj∈[0,10]cUj: the upper bound of investment cost of engineering characteristic j, cUj∈[0,10]pj: the required resource corresponding to engineering characteristic jB: the maximum total budgetP: the maximum total resources.


### 2.1. Objective Function and Decision Variables

The objective function in this model is to maximize the overall customer satisfactionCS, just as in formula (1). Whether the engineering characteristics are necessary can be viewed as 0 or 1 integer variables. If engineering characteristicj is selected, Nj is taken as 1. Otherwise, Nj is taken as 0. Another decision variable gj, the actual fulfillment level for engineering characteristic j, is directly influenced by cj, the corresponding investment cost. The detailed relationship expressions will be introduced in next section. Then these fulfillment levels of selected engineering characteristics will be synthesized according to their weights, to get the overall customer satisfaction. A traditional calculation method to get the weight of engineering characteristic j is as follows:(8)ρj′′=∑i=1mwi·Rij,where ρj′′ is the initial importance weight of the engineering characteristic j, wi is the relative importance weight of the requirements i, and Rij is the relationship between the customer requirement i and the engineering characteristic j, Rij∈{±1,±3,±5}.Moreover, the correlationship between different engineering characteristics is taken into consideration as in the following:(9)ρj=∑k=1nρj′′·Wjkwith Wjk denoting the correlation between customer requirement j and k, Wjk∈{±1,±3,±5}.


### 2.2. Relationship between Investment Cost and Actual Fulfillment Level

In reality, in formula (7), the relationship between investment cost and actual fulfillment levels is too complex to describe quantitatively. In order to simplify the problem, the relationship was generally viewed as linear in some past research [6, 32]. Here, the relationship is updated by a new expression, which can represent linear, concave, convex, or mixed relationship type.Considering the current known conditions firstly(i)
the independent variablecj ranges from 0 to 10;(ii)
the dependent variablegj ranges from 0 to 100%;(iii)
whencj is taken as 0, gj is taken as 0. When cj is taken as 10, gj is taken as 100%.A preliminary expression is suggested with the help of an auxiliary functionf as follows:(10)gjcj=fcj-f0f10-f0.Here, a tangent function is chosen as the auxiliary functionf. Moreover, in order to express more relationship types in the domain of definition with only one formula, a shape parameter u is introduced. The final expression is as follows:(11)gjcj=arctan⁡cj-uj+arctan⁡ujarctan⁡10-uj+arctan⁡uj.The shape parameteru varies with different engineering characteristics. Actually, u can be understood as a special investment cost value. At this point, the fulfillment level will get the fastest growing; that is, the slope of g(c) is maximum.Whenu goes to +∞ or -∞, g(c) is approximately linear, just as in Figure 2(a), which is an ideal situation.Figure 2
The relationship between the actual fulfillment level and the corresponding investment cost. (a)u=+∞ or -∞, (b) u≥10, (c) u≤0, and (d) 0<u<10.
(a) (b) (c) (d)Whenu≥10,  g(c) is a concave function, just as in Figure 2(b), which represents the case that the actual fulfillment level grows faster and faster with the increase of investment cost.Whenu≤0,  g(c) is a convex function, just as in Figure 2(c), which represents the case that the actual fulfillment level grows slower and slower with the increase of investment cost.When0<u<10,  g(c) is concave at the beginning and convex afterward, just as in Figure 2(d). The investment cost value at the inflection point is equal to u. The fulfillment level will increase sharply around this point.


### 2.3. Capability, Cost, and Resource Constraints

The fulfillment level for every engineering characteristic is limited by actual capability of the enterprise.gLj and gUj are corresponding to the lowest and the highest levels achieved by the enterprise. The actual fulfillment level gj should be within this range, just as in formula (2).Formula (3) and formula (4) are the cost constraints. The investment cost cj of engineering characteristic j should not surpass the upper bound cUj, and the total investment cost should not exceed the total budget limit B.Resource for every engineering characteristic is assumed to be known. The total resource should not exceed the resource limit for product developmentP as in formula (5).


## 2.1. Objective Function and Decision Variables

The objective function in this model is to maximize the overall customer satisfactionCS, just as in formula (1). Whether the engineering characteristics are necessary can be viewed as 0 or 1 integer variables. If engineering characteristicj is selected, Nj is taken as 1. Otherwise, Nj is taken as 0. Another decision variable gj, the actual fulfillment level for engineering characteristic j, is directly influenced by cj, the corresponding investment cost. The detailed relationship expressions will be introduced in next section. Then these fulfillment levels of selected engineering characteristics will be synthesized according to their weights, to get the overall customer satisfaction. A traditional calculation method to get the weight of engineering characteristic j is as follows:(8)ρj′′=∑i=1mwi·Rij,where ρj′′ is the initial importance weight of the engineering characteristic j, wi is the relative importance weight of the requirements i, and Rij is the relationship between the customer requirement i and the engineering characteristic j, Rij∈{±1,±3,±5}.Moreover, the correlationship between different engineering characteristics is taken into consideration as in the following:(9)ρj=∑k=1nρj′′·Wjkwith Wjk denoting the correlation between customer requirement j and k, Wjk∈{±1,±3,±5}.


## 2.2. Relationship between Investment Cost and Actual Fulfillment Level

In reality, in formula (7), the relationship between investment cost and actual fulfillment levels is too complex to describe quantitatively. In order to simplify the problem, the relationship was generally viewed as linear in some past research [6, 32]. Here, the relationship is updated by a new expression, which can represent linear, concave, convex, or mixed relationship type.Considering the current known conditions firstly(i)
the independent variablecj ranges from 0 to 10;(ii)
the dependent variablegj ranges from 0 to 100%;(iii)
whencj is taken as 0, gj is taken as 0. When cj is taken as 10, gj is taken as 100%.A preliminary expression is suggested with the help of an auxiliary functionf as follows:(10)gjcj=fcj-f0f10-f0.Here, a tangent function is chosen as the auxiliary functionf. Moreover, in order to express more relationship types in the domain of definition with only one formula, a shape parameter u is introduced. The final expression is as follows:(11)gjcj=arctan⁡cj-uj+arctan⁡ujarctan⁡10-uj+arctan⁡uj.The shape parameteru varies with different engineering characteristics. Actually, u can be understood as a special investment cost value. At this point, the fulfillment level will get the fastest growing; that is, the slope of g(c) is maximum.Whenu goes to +∞ or -∞, g(c) is approximately linear, just as in Figure 2(a), which is an ideal situation.Figure 2
The relationship between the actual fulfillment level and the corresponding investment cost. (a)u=+∞ or -∞, (b) u≥10, (c) u≤0, and (d) 0<u<10.
(a) (b) (c) (d)Whenu≥10,  g(c) is a concave function, just as in Figure 2(b), which represents the case that the actual fulfillment level grows faster and faster with the increase of investment cost.Whenu≤0,  g(c) is a convex function, just as in Figure 2(c), which represents the case that the actual fulfillment level grows slower and slower with the increase of investment cost.When0<u<10,  g(c) is concave at the beginning and convex afterward, just as in Figure 2(d). The investment cost value at the inflection point is equal to u. The fulfillment level will increase sharply around this point.


## 2.3. Capability, Cost, and Resource Constraints

The fulfillment level for every engineering characteristic is limited by actual capability of the enterprise.gLj and gUj are corresponding to the lowest and the highest levels achieved by the enterprise. The actual fulfillment level gj should be within this range, just as in formula (2).Formula (3) and formula (4) are the cost constraints. The investment cost cj of engineering characteristic j should not surpass the upper bound cUj, and the total investment cost should not exceed the total budget limit B.Resource for every engineering characteristic is assumed to be known. The total resource should not exceed the resource limit for product developmentP as in formula (5).


## 3. Solving Algorithm

The detailed G-ICA algorithm steps are shown in Figure3.Figure 3
G-ICA algorithm flow chart.Step 1.
A hybrid real-binary concatenated coding method is applied for the encoding and decoding of a solution.Step 2.
A penalty function method is used as the constraints handling method.Step 3 (empires construction).
The initial countries are generated and their corresponding objective function values are calculated. According to the power value of each country, which is a normalized expression of the objective function value, the most powerful countries are chosen as imperialists and the remained other countries are chosen as colonies. After the colonies are assigned to every imperialist randomly, the empires are constructed.Step 4 (assimilation).
Assimilation occurs inside every empire. According to the basic ICA, imperialists try to assimilate their colonies and make them similar to themselves, which aims to make solutions represented by colonies tending to solutions represented by imperialists. Here, the assimilation of real variables is the same as the original ICA algorithm.
But for integer variables, a certain modification has to be conducted. A multipoint crossover operator and a multipoint mutation operator are introduced, which refers to GA. The crossover probability isPcr, whose value is generally larger than 0.6. The mutation probability is Pmu, which is generally taken less than 0.1.Step 5 (revolution).
In each decade, revolution operators are performed on some of the colonies in each empire. There are many ways to realize revolution, such as swap, insertion, reversion, and perturbation. In this paper, some colonies are replaced by newly generated colonies. The population percentage that revolts is shown byPre.Step 6 (imperialist update).
After a series of operators inside every empire, the cost of every imperialist should be updated if there exists any colony whose cost is lower than the imperialist’s. In this case, the colony will become the imperialist in the current empire and vice versa.Step 7.
Calculate the total cost of each empire.Step 8 (empires competition).
The weakest colony of all empires will be selected and assigned to some empire. The lucky empire is chosen randomly; that is to say, the strongest empire does not always gain the initiative.Step 9 (elimination of the powerless empire).
When an empire loses all of the colonies, the empire will be eliminated from the population.Step 10 (stop criterion).
The stopping criterion is considered when there is only one empire that governs all of the colonies or the decades reached the limitations.


## 4. Case Studies and Discussions

### 4.1. Numerical Examples and Discussions

The G-ICA is used to solve 7 MINLP problems taken from different literatures. They are listed in the Appendix. The number of unknown decision variables in these problems varies from 2 to 10. More salient features of the test problems are listed in Table1.Table 1
Salient features of the test problems.
Problem ⁢Number of variables Objective function is nonlinear or not Number of constraints Total Real Integer Total Bounds Inequality Equality Nonlinear Linear Nonlinear Linear 1 2 1 1 Nonlinear 3 2 1 0 0 0 2 2 1 1 Linear 4 2 1 1 0 0 3 3 2 1 Nonlinear 6 3 1 2 0 0 4 3 2 1 Linear 4 2 0 1 1 0 5 4 2 2 Nonlinear 7 4 0 2 1 0 6 7 3 4 Nonlinear 16 7 4 5 0 0 7 10 5 5 Nonlinear 6 5 0 0 0 1Performance of G-ICA is compared with GA, (μ+λ)-ES, (μ,λ)-ES, and standard simulated annealing algorithm (SSA). Another two improved algorithms for MINLPs, RST2U, and RST2ANU algorithm, mentioned in Section 1, also take part in the comparison. All the algorithms were programmed by MATLAB R2014b, and these test problems are done on a personal computer equipped with Intel Core i7-5500U CUP, 2.4 GHz processor, 4 G RAM under Windows 8.1 operating system. The results obtained are given in Table 2.Table 2
Results obtained by using the following algorithm.
Problem  G-ICA GA (μ+λ)-ES (μ,λ)-ES SSA RST2U RST2ANU 1 p 100 90 100 100 100 100 100 t 0.0268 0.5294 0.2467 5.7895 0.0998 0.0325 0.0301 FEs 162 8360 1667 38200 1227 372 392  2 p 90 100 60 100 90 100 100 100 t 0.0413 0.0627 0.5872 0.4972 2.4193 0.1660 0.0767 0.0798 FEs 1016.2 1514.6 19915 2857 16250 2113 1132 1048  3 p 42 100 1 20 7 2 66 64 t 6.0803 18.0082 9.9711 36.9382 55.5372 29.325 26.8132 26.8426 FEs 303480 1000263 330090 339173 376063 319153 303552 303273  4 p 65 1 5 1 23 26 29 t 5.4450 70.4521 230.7150 266.5970 30.3212 19.0057 20.0926 FEs 260078 1992365 1852329 2281264 372594 253928 267426  5 p 59 — 51 — 39 35 39 t 55.4304 — 156.1188 — 453.7673 615.7494 629.7191 FEs 5000049 — 1061945 — 5135906 5247020 5132926  6 p 37 — 9 1 1 89 89 t 5.3581 — 144.6187 187.3324 38.2674 12.4025 12.5823 FEs 450103 — 807780 1007993 500240 100100 91705  7 p 24 — — — — 3 3 t 36.7097 — — — — 301.2364 297.8793 FEs 5000259 — — — — 5045722 5012452 Index:p—the percentage of success, t—average time in seconds used by the algorithm in achieving the optimal solution in the case of the successful runs, and FEs—the average number of function evaluations in the case of successful runs.Each problem is executed 100 runs with all the seven algorithms. Each run is initiated using a different set of initial population. A run is considered a success if achieved value of the objective function is within 1% of the known optimal value (in case the optimal value of the objective is zero, a run is considered success if the achieved absolute value of the objective function is less than 0.01). For each problem, the percentage of success (obtained as the ratio of the number of successful runs to total number of runs), the average number of function evaluations (FEs) in the case of successful runs, and the average time in seconds used by the algorithm in achieving the optimal solution in the case of the successful runs are also listed.From the results listed in Table2, it can be seen that, for most problems, the G-ICA can get an optimal solution in a much shorter time. For Problem 1, most algorithms can get 100% success rate, but G-ICA’s running time is the shortest and FEs is the smallest. For Problems 4, 5, and 7, on the premise of similar FEs, G-ICA can get a higher success rate in a much shorter time. For Problem 2, when the maximum number of function evaluations (MAXFEs) is 1000, G-ICA gets 90% success rate, which is lower than RST2U and RST2ANU. But when the MAXFEs reaches 1500, G-ICA gets 100% success rate, whose running time is still shorter than RST2U and RST2ANU. The same situation occurs in solving Problem 3.An interesting phenomenon is found in solving Problem6. The global optimal solution of Problem 6 given in the Appendix is x1,x2,x3,y1,y2,y3,y4,f=(0.2,1.280624,1.954483,1,0,0,1,3.557463). And the 37 successful runs of G-ICA are all around the global optimal solution. But RST2U and RST2ANU find another local optimal solution in searching process, x1,x2,x3,y1,y2,y3,y4,f=(0.2,0.8,1.9079,1,1,0,1,3.5796). In the 89 successful runs of RST2ANU, there are about 70% points around the local optimal solution, and the following are around the global optimal solution. In our calculating method of percentage of success, all these points are viewed as successful. So the data in Table 2 shows that G-ICA is far behind in the percentage of success index. This may be an evidence to prove that G-ICA can guarantee global optimal.In fact, the number of initial populations has a great influence on the performances of algorithms, which has been proposed by Engelbrecht [33]. In this paper, each run for each algorithm is initiated using a different set of initial population, and every index’s average value is taken as reference, in order to compare these algorithms fairly. And for G-ICA, both the number of initial countries Npop and number of imperialists Nimp, that is, the number of groups which the initial countries are divided into, are set differently. Here, an experiment is done for Problem 2 to research how these parameters influence the algorithm’s performances. The iteration stop criteria of G-ICA is to reach the maximum FEs. The results are shown in Figures 4 and 5.Figure 4
The influence ofNpop on percentage of success (a), average FEs of successful runs (b), and average time of the successful runs (c) for Problem 2.
(a) (b) (c)Figure 5
The influence ofNimp on percentage of success (a), average FEs of the successful runs (b), and average time of successful runs (c) for Problem 2.
(a) (b) (c)For the percentage of success, it can be seen from Figures4(a) and 5(a) that more FEs provide higher percentage of success. But the parameters Npop and Nimp’s influences to this index are hard to describe, because the relationship between Npop (Nimp) and the percentage of success is wavy.From the view of index FEs, the performance has become deteriorative with the increase ofNpop, just as in Figure 4(b), because the average FEs increase. But the average FEs will decrease when Nimp becomes larger, just as in Figure 5(b).The parameters’ influences to the average running time are shown clearly in the figures. On the premise of the same maximum FEs, the largerNpop is, the fewer iteration times are, so the average time is shorter, just as in Figure 4(c). But the larger Nimp is, the more times of empires competitions are, so the average time is longer, just as Figure 5(c). In conclusion, if you want to pursue a shorter running time, you may set a larger Npop and a smaller Nimp.


### 4.2. Case of the Washing Machine and Discussions

The data of this case is adapted from a washing machine case [13, 32, 34]. The QFD of a washing machine is as in Figure 6.Figure 6
QFD matrix for a washing machine.The importance weight matrix of engineering characteristics is as follows:(12)ω=54322T.The relationship matrix between requirements and engineering characteristics is as follows:(13)R=5015535110301513113550510.The correlationship matrix between engineering characteristics is as follows:(14)W=1033901000301333031190311.The initial importance weight of engineering characteristicj is calculated by formula (8). And the result is as follows:(15)ρ′′=6222245238T.The final importance weight of engineering characteristicj is calculated by formula (9). And the result is as follows:(16)ρ=W·ρ′′=63222480348720T.Here, the corresponding shape parameters are −600000, 200000,-1.3, and −0.5, 2.The object function of the optimization problem is as the following formula:(17)maxCS=632·N1·g1+22·N2·g2+480·N3·g3+348·N4·g4+720·N5·g5Subjectto0.1≤g1≤0.5,0.1≤g2≤0.4,0.2≤g3≤0.9,0.2≤g4≤0.9,0.2≤g5≤0.7,0≤c1≤2,0≤c2≤3,0≤c3≤6,0≤c4≤5,0≤c5≤4,N1·c1+N2·c2+N3·c3+N4·c4+N5·c5≤15,4N1+2N2+4N3+4N4+4N5≤17,g1=arctan⁡c1+600000+arctan⁡-600000arctan⁡10+600000+arctan⁡-600000,g2=arctan⁡c2-200000+arctan⁡200000arctan⁡10-200000+arctan⁡200000,g3=arctan⁡c3+1.3+arctan⁡-1.3arctan⁡10+1.3+arctan⁡-1.3,g4=arctan⁡c4+0.5+arctan⁡-0.5arctan⁡10+0.5+arctan⁡-0.5,g5=arctan⁡c5-2+arctan⁡2arctan⁡10-2+arctan⁡2.The model is solved using G-ICA and RST2ANU in the same computer and every algorithm runs 10 times. The results are as the first row in Table3. On the premise of getting similar solutions, the average FEs of RST2ANU are much smaller than that of G-ICA, but G-ICA’s average time is still shorter. All of these two algorithms’ solutions are clustered together. And one solution is taken as the decision reference: (18)c1,c2,c3,c4,c5,N1,N2,N3,N4,N5=1.99,2.31,5.54,4.53,2.80,1,0,1,1,1.Table 3
The results of washing machine case using G-ICA and RST2ANU.
n Number of variables G-ICA (s) RST2ANU (s) t 2 / t 1 Total Real Integer t 1(s) FEs1 Mean1 Std1 t 2(s) FEs2 Mean2 Std2 1 5 10 5 5 3.4 2 × 1 0 6 1373.6 0.8570 29.9 2 × 1 0 5 1373.1 0.7462 8.79 2 10 20 10 10 130.1 2 × 1 0 8 2714.9 12.8653 972.3 1.5 × 1 0 7 2716.4 4.5474 7.47 3 15 30 15 15 1180.4 2 × 1 0 10 4043.1 10.8875 9585.5 1.8 × 1 0 9 4044.2 11.2231 8.12 4 25 50 25 25 23609.23 2 × 1 0 13 8653.68 — 118050.15 2 × 1 0 12 8556.23 — 5.00 Index:n—the number of possible engineering characteristics, t—the average running time in seconds, FEs—the average number of function evaluations, mean—the average objective function value, and Std—the standard deviations of objective function values.The solution means that, considering the limited resource in this design phase, the engineering characteristic “noise level” can be ignored for the time being. The investment costs of the other four engineering characteristics are 1.99, 5.54, 4.53, and 2.80, respectively. And the overall customer satisfaction is 1375.4.Beyond that, another three cases adapted from the above case are also implemented. 10, 15, and 25 possible engineering characteristics are assumed in cases 2, 3, and 4. The other problem parameters are generated randomly. Cases 2 and 3 are implemented 10 times, but case 4 is only implemented 1 time because of its relatively long running time. The results are listed in Table3, and the time advantage of G-ICA is very obvious.


## 4.1. Numerical Examples and Discussions

The G-ICA is used to solve 7 MINLP problems taken from different literatures. They are listed in the Appendix. The number of unknown decision variables in these problems varies from 2 to 10. More salient features of the test problems are listed in Table1.Table 1
Salient features of the test problems.
Problem ⁢Number of variables Objective function is nonlinear or not Number of constraints Total Real Integer Total Bounds Inequality Equality Nonlinear Linear Nonlinear Linear 1 2 1 1 Nonlinear 3 2 1 0 0 0 2 2 1 1 Linear 4 2 1 1 0 0 3 3 2 1 Nonlinear 6 3 1 2 0 0 4 3 2 1 Linear 4 2 0 1 1 0 5 4 2 2 Nonlinear 7 4 0 2 1 0 6 7 3 4 Nonlinear 16 7 4 5 0 0 7 10 5 5 Nonlinear 6 5 0 0 0 1Performance of G-ICA is compared with GA, (μ+λ)-ES, (μ,λ)-ES, and standard simulated annealing algorithm (SSA). Another two improved algorithms for MINLPs, RST2U, and RST2ANU algorithm, mentioned in Section 1, also take part in the comparison. All the algorithms were programmed by MATLAB R2014b, and these test problems are done on a personal computer equipped with Intel Core i7-5500U CUP, 2.4 GHz processor, 4 G RAM under Windows 8.1 operating system. The results obtained are given in Table 2.Table 2
Results obtained by using the following algorithm.
Problem  G-ICA GA (μ+λ)-ES (μ,λ)-ES SSA RST2U RST2ANU 1 p 100 90 100 100 100 100 100 t 0.0268 0.5294 0.2467 5.7895 0.0998 0.0325 0.0301 FEs 162 8360 1667 38200 1227 372 392  2 p 90 100 60 100 90 100 100 100 t 0.0413 0.0627 0.5872 0.4972 2.4193 0.1660 0.0767 0.0798 FEs 1016.2 1514.6 19915 2857 16250 2113 1132 1048  3 p 42 100 1 20 7 2 66 64 t 6.0803 18.0082 9.9711 36.9382 55.5372 29.325 26.8132 26.8426 FEs 303480 1000263 330090 339173 376063 319153 303552 303273  4 p 65 1 5 1 23 26 29 t 5.4450 70.4521 230.7150 266.5970 30.3212 19.0057 20.0926 FEs 260078 1992365 1852329 2281264 372594 253928 267426  5 p 59 — 51 — 39 35 39 t 55.4304 — 156.1188 — 453.7673 615.7494 629.7191 FEs 5000049 — 1061945 — 5135906 5247020 5132926  6 p 37 — 9 1 1 89 89 t 5.3581 — 144.6187 187.3324 38.2674 12.4025 12.5823 FEs 450103 — 807780 1007993 500240 100100 91705  7 p 24 — — — — 3 3 t 36.7097 — — — — 301.2364 297.8793 FEs 5000259 — — — — 5045722 5012452 Index:p—the percentage of success, t—average time in seconds used by the algorithm in achieving the optimal solution in the case of the successful runs, and FEs—the average number of function evaluations in the case of successful runs.Each problem is executed 100 runs with all the seven algorithms. Each run is initiated using a different set of initial population. A run is considered a success if achieved value of the objective function is within 1% of the known optimal value (in case the optimal value of the objective is zero, a run is considered success if the achieved absolute value of the objective function is less than 0.01). For each problem, the percentage of success (obtained as the ratio of the number of successful runs to total number of runs), the average number of function evaluations (FEs) in the case of successful runs, and the average time in seconds used by the algorithm in achieving the optimal solution in the case of the successful runs are also listed.From the results listed in Table2, it can be seen that, for most problems, the G-ICA can get an optimal solution in a much shorter time. For Problem 1, most algorithms can get 100% success rate, but G-ICA’s running time is the shortest and FEs is the smallest. For Problems 4, 5, and 7, on the premise of similar FEs, G-ICA can get a higher success rate in a much shorter time. For Problem 2, when the maximum number of function evaluations (MAXFEs) is 1000, G-ICA gets 90% success rate, which is lower than RST2U and RST2ANU. But when the MAXFEs reaches 1500, G-ICA gets 100% success rate, whose running time is still shorter than RST2U and RST2ANU. The same situation occurs in solving Problem 3.An interesting phenomenon is found in solving Problem6. The global optimal solution of Problem 6 given in the Appendix is x1,x2,x3,y1,y2,y3,y4,f=(0.2,1.280624,1.954483,1,0,0,1,3.557463). And the 37 successful runs of G-ICA are all around the global optimal solution. But RST2U and RST2ANU find another local optimal solution in searching process, x1,x2,x3,y1,y2,y3,y4,f=(0.2,0.8,1.9079,1,1,0,1,3.5796). In the 89 successful runs of RST2ANU, there are about 70% points around the local optimal solution, and the following are around the global optimal solution. In our calculating method of percentage of success, all these points are viewed as successful. So the data in Table 2 shows that G-ICA is far behind in the percentage of success index. This may be an evidence to prove that G-ICA can guarantee global optimal.In fact, the number of initial populations has a great influence on the performances of algorithms, which has been proposed by Engelbrecht [33]. In this paper, each run for each algorithm is initiated using a different set of initial population, and every index’s average value is taken as reference, in order to compare these algorithms fairly. And for G-ICA, both the number of initial countries Npop and number of imperialists Nimp, that is, the number of groups which the initial countries are divided into, are set differently. Here, an experiment is done for Problem 2 to research how these parameters influence the algorithm’s performances. The iteration stop criteria of G-ICA is to reach the maximum FEs. The results are shown in Figures 4 and 5.Figure 4
The influence ofNpop on percentage of success (a), average FEs of successful runs (b), and average time of the successful runs (c) for Problem 2.
(a) (b) (c)Figure 5
The influence ofNimp on percentage of success (a), average FEs of the successful runs (b), and average time of successful runs (c) for Problem 2.
(a) (b) (c)For the percentage of success, it can be seen from Figures4(a) and 5(a) that more FEs provide higher percentage of success. But the parameters Npop and Nimp’s influences to this index are hard to describe, because the relationship between Npop (Nimp) and the percentage of success is wavy.From the view of index FEs, the performance has become deteriorative with the increase ofNpop, just as in Figure 4(b), because the average FEs increase. But the average FEs will decrease when Nimp becomes larger, just as in Figure 5(b).The parameters’ influences to the average running time are shown clearly in the figures. On the premise of the same maximum FEs, the largerNpop is, the fewer iteration times are, so the average time is shorter, just as in Figure 4(c). But the larger Nimp is, the more times of empires competitions are, so the average time is longer, just as Figure 5(c). In conclusion, if you want to pursue a shorter running time, you may set a larger Npop and a smaller Nimp.


## 4.2. Case of the Washing Machine and Discussions

The data of this case is adapted from a washing machine case [13, 32, 34]. The QFD of a washing machine is as in Figure 6.Figure 6
QFD matrix for a washing machine.The importance weight matrix of engineering characteristics is as follows:(12)ω=54322T.The relationship matrix between requirements and engineering characteristics is as follows:(13)R=5015535110301513113550510.The correlationship matrix between engineering characteristics is as follows:(14)W=1033901000301333031190311.The initial importance weight of engineering characteristicj is calculated by formula (8). And the result is as follows:(15)ρ′′=6222245238T.The final importance weight of engineering characteristicj is calculated by formula (9). And the result is as follows:(16)ρ=W·ρ′′=63222480348720T.Here, the corresponding shape parameters are −600000, 200000,-1.3, and −0.5, 2.The object function of the optimization problem is as the following formula:(17)maxCS=632·N1·g1+22·N2·g2+480·N3·g3+348·N4·g4+720·N5·g5Subjectto0.1≤g1≤0.5,0.1≤g2≤0.4,0.2≤g3≤0.9,0.2≤g4≤0.9,0.2≤g5≤0.7,0≤c1≤2,0≤c2≤3,0≤c3≤6,0≤c4≤5,0≤c5≤4,N1·c1+N2·c2+N3·c3+N4·c4+N5·c5≤15,4N1+2N2+4N3+4N4+4N5≤17,g1=arctan⁡c1+600000+arctan⁡-600000arctan⁡10+600000+arctan⁡-600000,g2=arctan⁡c2-200000+arctan⁡200000arctan⁡10-200000+arctan⁡200000,g3=arctan⁡c3+1.3+arctan⁡-1.3arctan⁡10+1.3+arctan⁡-1.3,g4=arctan⁡c4+0.5+arctan⁡-0.5arctan⁡10+0.5+arctan⁡-0.5,g5=arctan⁡c5-2+arctan⁡2arctan⁡10-2+arctan⁡2.The model is solved using G-ICA and RST2ANU in the same computer and every algorithm runs 10 times. The results are as the first row in Table3. On the premise of getting similar solutions, the average FEs of RST2ANU are much smaller than that of G-ICA, but G-ICA’s average time is still shorter. All of these two algorithms’ solutions are clustered together. And one solution is taken as the decision reference: (18)c1,c2,c3,c4,c5,N1,N2,N3,N4,N5=1.99,2.31,5.54,4.53,2.80,1,0,1,1,1.Table 3
The results of washing machine case using G-ICA and RST2ANU.
n Number of variables G-ICA (s) RST2ANU (s) t 2 / t 1 Total Real Integer t 1(s) FEs1 Mean1 Std1 t 2(s) FEs2 Mean2 Std2 1 5 10 5 5 3.4 2 × 1 0 6 1373.6 0.8570 29.9 2 × 1 0 5 1373.1 0.7462 8.79 2 10 20 10 10 130.1 2 × 1 0 8 2714.9 12.8653 972.3 1.5 × 1 0 7 2716.4 4.5474 7.47 3 15 30 15 15 1180.4 2 × 1 0 10 4043.1 10.8875 9585.5 1.8 × 1 0 9 4044.2 11.2231 8.12 4 25 50 25 25 23609.23 2 × 1 0 13 8653.68 — 118050.15 2 × 1 0 12 8556.23 — 5.00 Index:n—the number of possible engineering characteristics, t—the average running time in seconds, FEs—the average number of function evaluations, mean—the average objective function value, and Std—the standard deviations of objective function values.The solution means that, considering the limited resource in this design phase, the engineering characteristic “noise level” can be ignored for the time being. The investment costs of the other four engineering characteristics are 1.99, 5.54, 4.53, and 2.80, respectively. And the overall customer satisfaction is 1375.4.Beyond that, another three cases adapted from the above case are also implemented. 10, 15, and 25 possible engineering characteristics are assumed in cases 2, 3, and 4. The other problem parameters are generated randomly. Cases 2 and 3 are implemented 10 times, but case 4 is only implemented 1 time because of its relatively long running time. The results are listed in Table3, and the time advantage of G-ICA is very obvious.


## 5. Conclusions and Future Works

Determining the engineering characteristics to be implemented among so many possible selections and their corresponding actual fulfillment levels is an important QFD decision problem. But actual engineering problems’ increasing complexity may lead to a huge QFD matrix, which is time-consuming to analyze and tackle it. In order to improve QFD decision efficiency, an improved algorithm combining ICA and GA, called G-ICA, is proposed and applied to solve a MINLP model abstracted from the practical decision problem.In the process of defining problems, the relationship between investment cost and fulfillment levels of engineering characteristics is focused. This relationship is generally very complicated in practice and is often simplified to linear. Here, a new relationship expression is proposed to express linear, concave, convex, and mixed type with different values of shape parameteru. This shape parameter is assumed to be determined by the enterprises according to their experience, but with no more discussions in this paper. In the future work, more influence factors in practice will be taken into account while establishing QFD models, in order to accord with the actual product development situations better.The G-ICA algorithm is compared with GA, (μ+λ)-ES, (μ,λ)-ES, SSA, RST2U, and RST2ANU by calculating 7 MINLP problems in the same personal computer. Three indexes, percentage of success p, average time of successful runs t, and average FEs in the case of successful runs, are used to evaluate these algorithms. The results show that, comparing with other algorithms, the G-ICA algorithm has a satisfactory global optimization and runtime performance in solving MINLP problems. Finally, 4 cases of washing machine with different scales are solved by G-ICA and RST2ANU. For every case, both algorithms’ optimal results are clustered around the same point and can be as references for designers, to determine engineering characteristics and actual fulfillment levels. Meanwhile, the results show the G-ICA’s better runtime performance than RST2ANU again.



---

*Source: 2601561-2016-11-07.xml*