Genetic Algorithm for Mixed Integer Nonlinear Bilevel Programming and Applications in Product Family Design

Many leader-follower relationships exist in product family design engineering problems. We use bilevel programming (BLP) to reflect the leader-follower relationship and describe such problems. Product family design problems have unique characteristics; thus, mixed integer nonlinear BLP (MINLBLP), which has both continuous and discrete variables and multiple independent lower-level problems, is widely used in product family optimization. However, BLP is difficult in theory and is an NP-hard problem. Consequently, using traditional methods to solve such problems is difficult. Genetic algorithms (GAs) have great value in solving BLP problems, and many studies have designed GAs to solve BLP problems; however, such GAs are typically designed for special cases that do not involve MINLBLP with one or multiple followers. Therefore, we propose a bilevel GA to solve these particularMINLBLP problems, which are widely used in product family problems.We give numerical examples to demonstrate the effectiveness of the proposed algorithm. In addition, a reducer family case study is examined to demonstrate practical applications of the proposed BLGA.


Introduction
With the evolution of the mass customization paradigm, product family has played an increasingly important role in modern production and has garnered significant attention.Product family optimization design includes productdesign-related knowledge, the product family structure, and customization design based on the same product platform to meet customer needs.Many leader-follower relationships exist in product family optimization design problems, for example, between platform and customization design [1], between modular design and product family architecture [2], and between product family and supply design [3,4].Consequently, many researchers have used this model to reflect the optimization relationship in product family design.Du et al. [5] formulated a Stackelberg game-theoretic model for joint optimization of product family configuration and scaling design, wherein a bilevel decision structure reveals coupled decision-making between module configuration and parameter scaling.Kristianto et al. [6] demonstrated the potential of a two-stage, bilevel stochastic programming framework for tackling engineer-to-order product customization.
Du et al. [7] proposed a leader-follower joint optimization model of product family configuration and supply chain design, and Ji et al. [8] studied green design modularity using bilevel optimization.In addition, Jiao et al. [9] proposed an underpinning decision structure that distinguishes a de facto leader-follower model rather than a single-level, all-in-one optimization problem.
Due to the impact of Stackelberg game theory, which was proposed by Stackelberg when researching marketing economics, researchers have studied bilevel programming (BLP) since the 1970s.Bracken and Gill proposed BLP models in 1973 and 1977, respectively.Candler and Norton proposed a formal definition of BLP and multilevel programming in their technology reports [10].BLP models are NP-hard problems.Prominent BLP results primarily include the following: the kth best method for special linear cases [11], replacing a lower-level problem with Karush-Kuhn-Tucker conditions to convert a bilevel problem to a single-level problem [12], and converting the problem to a single level by constructing a penalty function based on the duality gap [13].The monographs of Bard [14] and Dempe [15] present systematic surveys of BLP models.In leader-follower joint optimization for product family design, 0-1 mixed integer nonlinear bilevel BLP, which has multiple nonlinear and nonconvex low-level models, is involved.Therefore, traditional methods (e.g., the K-T method for convex bilevel BLP) are limited to solving specific types of the BLP models, which restrict their application.
The BLP used in product family design has unique characteristics.In this type of programming, both continuous and discrete variables such as 0-1 variables (the most commonly used) are employed.In product family design problems, more than one lower-level model, which are nonlinear and nonconvex to model, are required.In addition, in this type of BLP, the leader's constraints often contain follower variables.Mixed integer nonlinear BLP (MINLBLP) contains the above characteristics when applied to product family optimization.Note that MINLBLP combines integer programming and BLP.However, the discreteness of decision variables and multiple followers make problems more complex.Some researchers prefer to use intelligent algorithms to solve this problem.BLP with multiple lower-level models is more difficult to solve than with a model consisting of only one follower.Therefore, developing an intelligent algorithm to solve the MINLBLP used in product family design problems has significant value.
A genetic algorithm (GA) is a method to search for an optimal solution by simulating biological evolution (survival of the fittest).GAs are popular intelligent algorithms that have seen increasingly wide utilization in many fields.GAs have many advantages such as convergence and robustness.Thus, GAs are effective in solving optimization problems.Using a GA to solve a bilevel problem reduces the limitations which traditional methods have, which has been extensively studied.Consequently, many GA monographs for BLP exist.In 1998, Liu designed a GA for multilevel programming with multiple followers, in which high-level models do not contain the lowlevel models' decision variables [16].In 2002, Oduguwa and Roy used a GA to solve a BLP scheme that was developed to encourage limited asymmetric cooperation between two players [17].Kuo and Han applied bilevel linear programming to a supply chain distribution problem and developed an efficient method based on a hybrid of a GA and particle swarm optimization [18], and Lin et al. presented a modified GA to solve a bilevel linear programming network design problem [19].For a linear bilevel programming problem with interval coefficients in the upper-level objective, Fan and Li propose a genetic algorithm to obtain the best optimal solution and the worst one simultaneously [20].Li et al. convert an integer linear bilevel programming problem into a single-level binary implicit program and then develop an orthogonal genetic algorithm for solving such problem based on the orthogonal design [21].Based on a novel coding scheme, Li proposes a genetic algorithm with global convergence to solve nonlinear bilevel programming problems where the follower is a linear fractional program [22].However, even though GAs have been successfully applied to BLP, it is difficult for such GAs to solve MINLBLP in product family design problems.Such GAs focus on solving specific models for other problems.They can solve linear BLP, single-level programming, and BLP with multiple followers, in which leader constraints do not contain follower decision variables.However, these GAs are not suitable for solving MINLBLP, which contains 0-1 decision variables and multiple low-level models.A bilevel GA (BLGA) for MINLBLP is proposed to solve product family design optimization problems.
The proposed BLGA can handle MINLBLP with both continuous and discrete variables with one or multiple independent nonlinear and nonconvex lower-level models.The remainder of this paper is organized as follows.In Section 2, we pose the primary research question and derive MINLBLP.In Section 3, we propose the BLGA and describe its procedures.Numerical examples are provided in Section 3 to illustrate the proposed algorithm.In Section 4, we provide a case study of product family design to demonstrate an application of the proposed method.Conclusions are given in Section 5.

MINLBLP with One or Multiple
Independent Followers An MINLBLP model probably has one or multiple followers.
Based on the relationships among the multiple followers, such MINLBLP problems can be classified as dependent or independent and, for each particular case, there exists a specific model to describe it [23][24][25].The coupling relationships among the followers are not taken into account for most of the bilevel problems in product family design [26,27].Therefore, this paper focuses on the MINLBLP with independent followers.In an MINLBLP model with multiple independent followers, the followers are independent and have no coupling relationship, which means that the objective function and the set of constraints of each follower only include the leader's variables and its own variables [24].
Here, assume the leader first chooses its decision variables  = { where  means the number of followers.

BLGA Design
3.1.Assumption.Analysis of feasibility of the proposed algorithm illustrates that two assumptions must be satisfied to obtain good performance.These assumptions ensure the effectiveness of the proposed BLGA by restraining the algorithm's process.
The first assumption is that solutions must meet the requirements of BLP, which means that the solutions must satisfy the constraints of both the upper-and lower-level programs.Meaningful solutions can only be obtained if the algorithm is feasible.
We propose a BLGA that meets these constrain requirements.First, the binary code encodes the leader variables and initializes them according to their bound constraints.Second, by taking the initial leader variables as parameters, each follower programming uses the lower-level GA to obtain the best solution in the last generation.Note that the best solution in the last generation provided by the lower-level GA is a good approximation for the optimal solution, which meets the lower-level constraints.Then the lower-level approximation solutions are brought back to the upper-level GA to solve the leader program.In the upper-level GA, the solutions should be verified whether they satisfy the leader's constraints.
The second assumption is that the algorithm has convergence.As the proposed BLGA meets the first condition, it can obtain the optimal solution only if it has convergence, which is an important measurement of computing capacity.Convergence of a GA means that the last cycle of a finite series of cycles can yield the optimal solution.This assumption ensures that, after comparing feasible solutions, the optimal solution can be found.An algorithm is only effective and feasible when it satisfies this condition.
Currently, many GAs have been designed to solve special types of problems that have convergence [28][29][30].Markov chains are often used to evaluate convergence.If a GA whose mutation rate is fixed and greater than 0 takes an elitist strategy and the binary code strategy, it has convergence that cannot be influenced by the initial population.Using the elitist strategy allows the best chromosome to be removed from operations in order to avoid disappearing in the operations.Moreover, many factors influence convergence, such as the maximal number of generations, crossover rate, mutation rate, and crossover and mutation methods.
According to the above assumptions, even though the best solution provided by the BLGA does not guarantee to be the optimal solution, it may be a good approximation for the global optimal solution in reasonable time.We can calculate the problem several times, using different seeds, to find the closest approximation solutions.Numerical examples are given in Section 3.4 to verify the algorithm.

Design Process of BLGA.
On the basis of the above assumptions, this study proposes a BLGA to solve MINLBLP with one or multiple independent followers that can improve efficiency in generating solutions.The essential process is as follows.First, leader decision variables are initialized according to their bound constraints.Second, each lower-level GA takes leader decision variables  = { 1 ,  2 , . . .,   , . . .,   } as parameters to solve the th follower problem to determine the best solution  *  () and the value of  *  (,  *  ()), where  = 1, . . ., .Third, the best solution  *  () for each follower and the leader decision variables  is returned to the leader problem to determine if they satisfy other leader constraints.Their fitness values are then evaluated.Fourth, selection, crossover, and mutation operations are performed for these leader decision variables.Fifth, this cycle continues until maximal iterations are reached.The best optimal solutions  * and  * 1 , . . .,  *  are then recorded.However, there is no guarantee that the optimal solutions will be the best solution for MINLBLP with one or multiple independent followers because the GA can only ensure that solutions belong to the constraint region rather than the inducible region.It is necessary to run this process several times to obtain several optimal solutions.These optimal solutions are then compared to determine the best solution.
The process of the proposed BLGA is shown in Figure 1.

BLGA Procedure.
The procedure of the proposed BLGA is as follows.
Step 1 (verify if the process is terminated).The termination condition is that the number of runs has reached the maximal number.To determine if the termination condition has been satisfied, the maximal number of runs  must be input.If this condition is satisfied, each optimal solution is compared to determine the best solution.If the termination condition is not satisfied, the process proceeds to Step 2.
Step 2 (initialize the leader population with leader population size   ).According to their bound constraints, the leader population is initialized and these variables are encoded.Two strategies are used for the encoding of variables.The proposed method uses the binary code strategy, which uses binary vectors as chromosomes, to represent the values of continuous variables with unequal lengths.The variable encoding length  is expressed as follows: where  max and  min represent the upper and lower bounds of the variables, respectively, and  is precision.The proposed method also uses the real code strategy for encoding, which uses ordered arrays to represent discrete variable values.
Step 3 (run lower-level GAs).Each follower programming takes the initial leader population  = { Step 5 (verify if the cycle is terminated).The termination condition is that the number of generations has reached the maximal number.If the maximal number is reached, the chromosome with the greatest fitness value is considered the best solution.The best leader chromosome  * and each best follower chromosome  *  (),  = 1, . . ., , for each generation are recorded and the process returns to Step 1.If the maximal number has not been reached, the process continues to Step 6.Note that inputting a correct maximal number of generations is very important to an algorithm's convergence.If the number is too small, the algorithm may not find the best solution after termination of the cycle.In contrast, if the number is too large, efficiency may be reduced.A careful and deliberate choice is required.
Step 6 (selection operation).The roulette wheel selection method is used to select   −1 good chromosomes as parents for the crossover operation.Simultaneously, another good chromosome is selected as the best chromosome within the current population according to the elitist strategy.
Step 7 (crossover operation).Input the crossover rate   .Each of the parents selected in Step 6 is picked randomly, and each pair uses a uniform crossover method to undergo crossover with the given probability.Note that the new population includes these new chromosomes.
Step 8 (mutation operation).The proposed method uses two mutation methods.For continuous variables that use the binary code strategy for encoding, the mutation rate   is input and the uniform mutation method with the given probability is used.For discrete variables that use the real code strategy, new chromosomes are formulated as follows: where   and  ∧  are the individuals before and after mutating, respectively, and lb  and ub  are the lower and upper bounds of the integer variables, respectively.After completion of Step 8, the process returns to Step 3.

Numerical Examples.
The proposed BLGA for solving MINLBLP with one or multiple independent followers is realized using MATLAB.To verify the effectiveness and robustness of the proposed BLGA, we give numerical examples obtained with a personal computer with the following parameters: population size is 50; maximal number of generations is 200; precision of binary code is 0.01; crossover rate is 0.8; mutation rate is 0.01; number of experiments is fifteen.The results were compared with the exact solutions obtained from traditional methods for solving MINLBLP to test the quality of the computed solution.The best, worst, average, and median results and standard deviation value are presented.
Example 1 (leader and follower decision variables are discrete).
Example 1 is a simple integer nonlinear BLP problem that can be solved using the enumeration method.The best, worst, average, and median results and standard deviation value are shown in Table 1.
Table 1 shows the results of Example 1.The standard deviation value of  * ( * ,  * ) is zero, which means that each run of this BLGA in this test finds the same upper objective function value.However, the standard deviation of  * ( * ,  * ) is 3.04, which means that there exists more than one optimal solution found by the proposed BLGA in this example.From the result of each run, we find that the solutions  = (2, 0) and  = (1, 1) are the optimal solutions.However, even though both of these solutions have the same upper optimal solution, their lower optimal solutions are different.When the upper-level solution is  = (2, 0), the lower-level result is −2.When the upper-level solution is  = (1, 1), the lower-level result is −8.The upper solution  = (1, 1) provides a better value for the lower-level objective function.Therefore, the best solution for this test is  = (1, 1) and  = (0, 2).The iterative process is shown in Figure 2. As seen in Figure 2(a), the maximal number of generations does not need to be 200 because the optimal solution has been determined in the first generation.Example 2 is simple continuous BLP that can be solved by the graphic method.The best, worst, average, and median results and standard deviation value are shown in Table 2.
As can be seen, the difference between the values of upper-level objective functions of the best result calculated by the proposed BLGA and the exact method is 0.0030, which is a negligible difference.This demonstrates the effectiveness of the proposed BLGA.In these 15 independent runs, except the worst result, the other results are the same.The iterative process of the best result is shown in Figure 3(a).As can be seen, the maximal number of generations does not need to be 200 because the optimal solution has been determined in the fifth generation.
The results calculated by the proposed BLGA and the exact method are shown in Table 3.
Table 3 shows the best and the worst results among these 15 runs, whose objective function values are the same as that calculated by the exact method, even though the solutions are different.The standard deviation value is near zero, which means that each independent run draws the similar objective function values.Figure 4(a) shows the iterative process of a best run, and the evolution processes for the upper-level and lower-level GAs are shown in Figure 4(b).
The results are shown in Table 4.In these 15 runs, there are 6 times in which the proposed BLGAs draw the different solutions with the best result.However, these differences are so small that can be negligible.Moreover, the best result is rather similar to the exact solution, which demonstrates the effectiveness of the proposed BLGA.The iterative process is shown in Figure 5(a), and the BLGA evolution processes for the upper-level and lower-level optimizations are shown in Figure 5(b).
The results are shown in Table 5.Table 5 shows the basic information about these 15 independent runs.Note that the differences between the best objective function values calculated by the proposed BLGA and the exact method are negligible.The iterative process is shown in Figure 6(a), and the BLGA evolution processes for the upper-level and lower-level optimizations are shown in Figure 6(b).

Algorithm Evaluation.
Based on the characteristic of MINLBLP with one or multiple independent followers in product family design (PFD), we have proposed a BLGA to solved this problem.We use five numerical examples to demonstrate the effectiveness and robustness of the proposed BLGA.We summarize the advantages of the proposed BLGA as follows.
First, MINLBLP with one or multiple independent followers has generality, which can illustrate many practical management and engineering problems, such as PFD.Many studies have examined using GAs to solve special BLP problems.Some algorithms require that leader constraints include no follower decision variables, which cannot solve the target MINLBLP problems.In addition, some algorithms require only continuous decision variables, which limit practical application.In theory, the proposed BLGA has no limitations on constraints and variables.
Second, the optimal solutions calculated by the proposed algorithm belong to the constraint region.Note that to ensure the feasibility of the optimal solutions, we have made the following assumptions prior to designing the proposed BLGA.In the BLGA, the fitness values of individuals that do not satisfy constraints are set to zero to avoid being selected as parents.This ensures that the optimal solutions satisfy all constrains, which belong to the constraint region.
Third, the proposed BLGA has good effectiveness and strong robustness, which has been demonstrated through the five numerical examples.Note that we use the binary code strategy to encode so that we can set the precision according to practical requirements.This strategy suits the mechanism of GAs, which can search for optimal solutions carefully.

Problem Description.
A case study of a simplified reducer family design is shown to demonstrate the proposed BLGA.A single-stage gear reducer includes three types of customized reducers, whose transmission ratios are 2, 3, and 4, respectively.To reduce cost and produce the three reducers simultaneously, we design platform parts and different customization parts, respectively.In this case, platform parts include a shell, a high-speed shaft, and a low-speed shaft, whose critical parameters are taken as the platform variables.The different customization parts include the high-speed gear and the low-speed gear, whose critical parameters are taken as the customization variables.The critical parameter of platform parts is shell size, which is described by the center distance and diameter of high-speed shaft   .The center distance is proportional to the total number of gears when the module is determined; thus, we take the total number of gears   as a platform variable.To optimize customization parts, we must optimize the selection of materials of gears   (a 0-1 variable) and the number of teeth of high-speed gear    ( = 1, 2, 3; the total number of gears is the sum of the number of lowspeed gear and high-speed gear teeth).Here,   = 1 indicates that the th reducer uses the jth material.
In the design process of the reducer family, we cannot consider platform optimization and customization optimization simultaneously.The platform designer who optimizes the platform aims to leverage enterprise profitability, which includes maximizing the performance of all reducers and minimizing total cost.The customization designers design customization variables to maximize the performance of each reducer.Note that these designers have different goals and decision variables; however, their decisions affect each other.This is a joint optimization problem.Compared with the customization designers, platform designers act as leaders because the platform is the basis of a product family, which is more important.The customization designers act as followers.Consequently, this case can be considered a Stackelberg game-based decision-making problem.

Model Development.
In this reducer family design problem, the leader is the platform design, whose decision variables are the total number of gears   and the diameter of high-speed shaft   .The objective is to maximize the performance of all reducers and minimize total cost.The constraints are the ranges of the decision variables.
In this problem, the performance of the reducer family is the sum of the performance of each reducer, which is influenced by bending resistance and abrasion resistance.We use the bending fatigue stress of high-speed gear    to define the bending resistance of the th reducer because bending resistance affects the bending and fracture of a gear and, compared with a low-speed gear, the high-speed gear is easier to break.Because it is influenced by contact fatigue stress, the abrasion resistance is defined by the contact fatigue stress    .In this simplified case, we describe the evolution of fatigue stress of the th reducer by the product of two types of fatigue stress.The evolution of fatigue stress of the reducer family is ∑       .This has an inverse relationship with the performance of the reducer; thus, the goal of maximizing the performance of all reducers can be replaced by minimizing the evolution of fatigue stress of the reducer family.
Here, the total cost is the sum of the costs of each reducer, which is influenced by manufacturing materials.Assume that different types of reducers that use the same type of material have the same production cost.The customization designer can select the manufacturing materials from metals A, B, and C. The production costs   ( = 1, 2, 3) are 50, 70, and 100, respectively.We define   as 0-1 variables that show whether the th reducer uses the jth material.The production cost of the th reducer is ∑     , and the production cost of the reducer family is ∑ ∑     .
We can describe the leader's objective of the upper-level model as minimizing the product between the evolution of fatigue stress and the total cost, which is given as follows: The lower-level models attempt to optimize the three types of reducers.For the th reducer, the decision variables are the diameter of low-speed shaft    , the number of teeth of high-speed gear    , and the materials selection variable   (,  = 1, 2, 3).The objectives of the lower-level models are to minimize the cost of reducer     .The constraints are the restrictions of fatigue stress, dimensions, ranges of indicators, ranges of variables, and material selections.
The constraints of fatigue stress are that bending fatigue stress and the contact fatigue stress are less than allowable stress, respectively, which are given as follows: where loading coefficient  is 1.4, region coefficient   is 2.5, tooth width   is 80 mm, and the torque of the high-speed shaft is 1.27 × 10 5 N⋅mm.The different material coefficients of metal are given in Table 6.The dimension constraints dictate that the diameters of the reference circles of low-speed and high-speed gears are greater than 20 mm longer than the diameters of the low-speed shaft and high-speed shaft, respectively, which are given as follows: Suppose the margin of error for the actual transmission ratio of the reducers is 10% above or below the nominal transmission ratio.The constraints of the indicator ranges include the ranges of the actual transmission ratio, which is expressed as follows: Each customization designer can only select a single material.This constraint is expressed as follows: In addition, there are also constraints for the ranges of variables.These constraints are given in Table 7.Note that The material coefficient of metal and the ranges of variables and indicators are shown in Tables 6 and 7, respectively.

Solution and Analysis of Results
. This model is an MINLBLP problem with three independent followers.We use the proposed BLGA to solve this problem with the following parameters: population size is 50; maximal number of generations is 300; precision of binary code is 0.01; and crossover rate is 0.8. Figure 7(a) shows that the best optimal solutions are found in the fiftieth generation.Figure 7(b) shows the evolution processes of the overall evaluation and cost.Figure 7(c) shows the overall evolution processes of the fatigue stress of the family and each reducer.The optimal results are given in Table 8, which shows the design of each reducer.

Conclusion and Further Study
MINLBLP can accurately model product family design, which has one or multiple nonlinear and nonconvex independent low-level models.Based on the characteristics of the model in product family optimization, the primary purpose of this study has been to propose a bilevel BLGA to solve product family optimization models.The proposed BLGA applies a GA to both the leader and followers and uses appropriate encoding strategies to handle continuous and discrete decision variables, which are common in product family optimization.The proposed BLGA overcomes many limitations of previously proposed GAs.For example, such GAs cannot solve a model in which leader constraints do not contain follower decision variables, which is common in product family optimization.The proposed BLGA overcomes this limitation.This paper models the joint optimization between product platform and product family designs by an MINLBLP approach in detail and presents a case study of a simplified reducer family design to demonstrate the feasibility and effectiveness of this model and BLGA.Currently, BLP is widely used in product family design optimization, which is a hot topic in engineering.GAs can solve BLP problems effectively.This paper only presents a specific application of a simplified reducer family design which describes the joint optimization between product platform design and product family design.In the future, we can do more research on other leader-follower relationships in the product family design problems and improve the BLGA efficiency by adjusting the encoding strategies and operation strategies.

Figure 3 (
b) shows the BLGA evolution processes for the upper-level and lower-level optimizations.Example 3 (leader decision variables are discrete and follower decision variables are continuous).min   (, ) = −
1 ,  2 , . . .,   , . . .,   } as a parameter.Note that the processes for each lower-level GA are the same.First, the follower population size   is input, and the follower population is initialized according to their bound constraints.The two strategies described in Step 2 are used to encode the population.Second, verify that the follower chromosomes satisfy the follower constraints.If the follower chromosome does not satisfy the requirement, set the fitness values of it to zero.The others are evaluating the fitness value based on the follower objective function.Third, the maximal number of generations is input to verify whether the current number of generations is equivalent to the maximal number.If the current generation reaches the maximal number, the best follower chromosome  *  (),  = 1, . . ., , is recorded and the solution's feasibility is verified.If the solution meets the lower-level constraints, it is returned to the leader program or else null is returned.If the current population does not reach the maximal number, selection, crossover, and mutation operations are performed and the process continues. = 1, . . ., , are null values.If at least one is null, the fitness value is set to zero.If no null values are found, (  ) are sorted in descending order to obtain the largest value, that is,  max .The current  max is compared with the previous  max , and the larger value is recorded as the new  max .The fitness function  max − ( * ,  * 1 , . . .,  *  ), which can ensure nonnegative fitness values, is used to evaluate the individuals and record their fitness values.If the constraints are not satisfied, the fitness value is set to zero.

Table 6 :
Material coefficients of metals.