Multi-Objective Optimization of Activation Time and Discharge Time of Thermal Battery Using a Genetic Algorithm Approach

Activation time and discharge time are important criteria for the performance of thermal batteries. In this work a heat transfer analysis is carried out on the working process of thermal batteries. The effects of the thicknesses of heat pellets which are divided into three groups and that of the thickness of insulation layers on activation time and discharge time of thermal batteries are numerically studied using Fluent 15.0 when the sum of the thickness of heating plates and insulation layers remain unchanged. According to the numerical results, the optimal geometric parameters are obtained by using multi-objective genetic algorithm. The results show that the activation time is mainly determined by the thickness of the bottom heat pellet, while the discharge time is determined by the thickness of the heat pellets and that of the insulation layers. The discharge time of the optimized thermal battery is increased by 4.08%, and the activation time is increased by 1.23%.


Introduction
Thermal batteries, which are widely used in missiles and other modern weapons, have the advantages of long storage time, no self-discharge, quick activation and high reliability [1][2][3][4]. Unlike conventional batteries, thermal batteries need to operate at very high temperatures (usually more than 350 • C [1,5]) in order to melt the electrolyte, which is a solid at room temperature, which causes them to activate and discharge outward. How to shorten the activation time and prolong the discharge time is a popular research topic that has attracted wide attention. Some researchers have studied the activation stage of thermal batteries, such as Kang et al. [6], who numerically simulated the activation of thermal batteries without a centre hole, and concluded that heat conduction is the main factor affecting the activation time. Witt et al. [7] used a 3-dimensional scanning laser Doppler vibrometer (3-D SLDV) to capture the temperature variation of thermal batteries during activation process. Their results indicated that the electrical properties of thermal batteries vary with the temperature. Koyuncu et al. [8] established a two-dimensional finite element model based on COMSOL Multiphysics to study the activation time of thermal batteries. The model provided a reliable heat design tool for the activation phase of thermal batteries. For the discharge stage, Krieger et al. [9,10] investigated the effects of gases and gas mixtures on global thermal conductivity values of porous thermal insulation packages of fully assembled low cost competent munition (LCCM) thermal batteries using transient and quasi-steady-state heat transfer techniques. Their results indicated that filling the thermal battery with an inert gas, such as krypton or xenon, significantly prolonged the discharge time of the thermal battery. Yang et al. [11] presented a one-dimensional mathematical model for a high 2 of 17 temperature lithium-aluminum, iron disulfide molten salt battery system. The model predicts that increasing the battery's temperature can reduce the precipitation of KCl and prolong the discharge time. Haimovich et al. [12] proposed the method of using an additional buffer salt to prolong the discharge time. Chen et al. [13] developed a thermal battery which greatly prolongs the discharge time compared to the traditional thermal battery with a low melting point electrolyte. Jeong et al. [14] conducted thermal analysis on large-capacity thermal batteries, reporting that adding additional heat pellets to the insulation layer lead to longer discharge time.
Genetic algorithms (GAS) are a class of stochastic algorithms well-suited for large-scale optimization, which have been successfully applied to optimize lithium batteries heating strategy [15], multi-physical field modeling [16], and battery shell structure design [17]. The rationale of the GA procedure is used to simulate the process of natural selection and mutation and determine the optimal solution to minimize the errors between simulated data and experimental data.
Although many studie have been conducted in the simulation of the discharge phase and activation phase of thermal batteries, most of these studies have not studied these two phases as a whole. In addition, there is no research on optimizing the activation time and discharge time of thermal battery. However, rapid response devices require thermal batteries to have a short activation time and have a long work time the thermal battery work. In general, short activation time and long discharge time are two essential characteristics of high-performance thermal batteries. In the manufacture of thermal battery, the thickness of heat pellets should be increased to shorten the activation time, while the thickness of insulators should be increased to extend the discharge time. Obviously, both cannot be satisfied at the same time for a given shell. Therefore, in the actual design of thermal battery, the activation stage and discharge stage should be considered comprehensively, and the optimized structure parameters should be gotten in order to achieve the best performance of thermal batteries. It is urgent to investigate the relationship between the thickness of heat pellets and insulators in order to obtain the best activation time and discharge time. This work uses ANSYS Fluent to analyze the activation process and discharge process of the thermal battery and compares the data with that given in the literature to verify the correctness of the simulation results. From the perspective of shortening the activation time and increasing the discharge time, a multi-objective optimization method based on genetic algorithm is proposed. The thickness of the heat pellets in the electroactive battery and that of the insulation layers on the upper and lower sides are used as the decision variables to optimize the activation time and the discharge time. Finally, the optimization results are verified by simulation.

Structural Composition and Parameters of Thermal Battery
In this paper we use the thermal battery simulation model reported by Jeong et al. [14], as shown in Figure 1a. The 13-stack thermal battery is composed of a case, thermal insulation, top and bottom stack regions, and the electroactive cells, which consist of 13 layers of unit cells and 14 layers of heat pellets (Fe/KClO 4 ). As shown in Figure 1c, the unit cell consists of cathode (FeS 2 ), electrolyte (LiF-LiCl-LiBr), anode (LiSi), and current collector (stainless steel). Stack regions with heat pellets and insulators A (mica), as shown in Figure 1b,d, are added to both ends of the electroactive cell to provide sufficient heat to both ends of the electroactive cell.

Governing Equations and Boundary Conditions
Assuming that the fluid is incompressible and the physical parameters are constant, the governing equations of the fluid are expressed as follows [18]: Continuity equation: Momentum equation: Energy equation: The boundary conditions are set as follows: Initial temperature of the battery: At the wall of the battery:

Governing Equations and Boundary Conditions
Assuming that the fluid is incompressible and the physical parameters are constant, the governing equations of the fluid are expressed as follows [18]: Continuity equation: Energy equation: The boundary conditions are set as follows: Initial temperature of the battery: At the wall of the battery: where h = 10 W m −2 K −1 [19], T am = 300 K, ε = 0.11 [19], σ = 5.67 W m −2 K −4 , n is the direction of outward normal. At the axis: Energies 2020, 13, 6477 4 of 17 The thermophysical parameters for all materials are listed in Table 1. The commercial software FLUENT 15.0 is used as the simulation tool in this study. The heat paper is defined as the air, and the other computational domains are defined as the solid domain. The time step size can be divided by three. The time step of 0 s to 4 s is 0.005 s, the time step of 4 s to 10 s is set as 0.1 s working as bumper to avoid the rapid change of time step size and after 10 s, the time step is set as 5 s, and the calculated time lasts until the end of the battery discharge time.

Setting of Heat Source
As previously mentioned, heat pellets burn to release energy, which melts the electrolytes when the thermal battery begins to discharge outwards. Assuming that heat pellets burn at the same time and produce heat evenly, the burning time ∆t hp = 0.42 s and the power density of a pellet q hp = 9.7957 × 10 9 w m −3 . When the electrolyte reaches its melting point (709.15 K) indicating that the thermal battery is activated. After the activation process, the Joule heating (exothermic) and chemical reaction heat (endothermic) are generated in the discharge phase. Joule heat is mainly generated at the cathode [14] and can be calculated by the following expression [21]: where R is the resistivity and i A is the current per cross sectional area (i.e., current density). The resistivity of each component is 1.8 × 10 −4 Ω m for cathode [21], and the current density is 500 A m −2 [20]. On the other hand, chemical reaction heat is mainly generated at the cathode and anode [14] and can be estimated as [21]: where i V is the current per unit volume (i.e., i V = i A /δ with δ being the radius of cell), T is the temperature, and ∂U 0 /∂T is the variation of the open-circuit cell potential with respect to temperature,

Grid Independence Tests
The thermal battery model is built by using Workbench 15.0, and the mesh is defined as a triangular grid. We test the mesh independence with four mesh systems of 229,379; 272,758; 317,859 and 811,353. Figure 2 shows the average temperature change curve of the first layer electrolyte with different mesh systems. When the mesh number increases, the average temperature T e1 of the first layer electrolyte gradually increases. When the mesh number exceeds 317,859, the value of T e1 tends to be consistent, and the maximum error is less than 0.1%. Therefore, 317,859 is selected as the number of grid division in the simulation calculation when the accuracy meets the requirements.
Energies 2020, 13, 6477 5 of 17 triangular grid. We test the mesh independence with four mesh systems of 229,379; 272,758; 317,859 and 811,353. Figure 2 shows the average temperature change curve of the first layer electrolyte with different mesh systems. When the mesh number increases, the average temperature e1 T of the first layer electrolyte gradually increases. When the mesh number exceeds 317,859, the value of e1 T tends to be consistent, and the maximum error is less than 0.1%. Therefore, 317,859 is selected as the number of grid division in the simulation calculation when the accuracy meets the requirements.

Correctness Verification of Numerical Simulation
In order to verify the correctness of the numerical simulation, we compare the calculated results with those of Jeong et al. [14]. Figure 3 shows the temperature curve of points 1 and 2 (as shown in Figure 1). It can be seen from Figure 3 that in the activation stage (0-1 s), the temperature in this paper is consistent with that in [14], so the activation time is not affected. In the discharge stage (after 1 s), the temperature in this paper was slightly higher than that reported in [14], which would slightly increase the discharge time obtained in this paper. Although there are differences between the simulation results and [14], the variation trend and error of the results (the maximum relative error is 9.64%) are within an acceptable range, and the optimization is carried out on the basis of the model established in this paper, that is to say, the optimization results are compared with the simulation results in this paper instead of [14]. Therefore, we conclude that the simulation error in this paper has no influence on the optimization result.

Correctness Verification of Numerical Simulation
In order to verify the correctness of the numerical simulation, we compare the calculated results with those of Jeong et al. [14]. Figure 3 shows the temperature curve of points 1 and 2 (as shown in Figure 1). It can be seen from Figure 3 that in the activation stage (0-1 s), the temperature in this paper is consistent with that in [14], so the activation time is not affected. In the discharge stage (after 1 s), the temperature in this paper was slightly higher than that reported in [14], which would slightly increase the discharge time obtained in this paper. Although there are differences between the simulation results and [14], the variation trend and error of the results (the maximum relative error is 9.64%) are within an acceptable range, and the optimization is carried out on the basis of the model established in this paper, that is to say, the optimization results are compared with the simulation results in this paper instead of [14]. Therefore, we conclude that the simulation error in this paper has no influence on the optimization result.

Effects of Heat Pellet Thickness and Insulation Layer Thickness on the Activation Time and Discharge Time
As shown in Figure 4, the heat pellets are divided into three groups. The 1st to 4th layers are defined as group a, the 5th to 10th layers are deemed as group b, and the 11th to 14th layers are considered as group c. It is stipulated that the thickness of heat pellets in the same group is the same,

Effects of Heat Pellet Thickness and Insulation Layer Thickness on the Activation Time and Discharge Time
As shown in Figure 4, the heat pellets are divided into three groups. The 1st to 4th layers are defined as group a, the 5th to 10th layers are deemed as group b, and the 11th to 14th layers are considered as group c. It is stipulated that the thickness of heat pellets in the same group is the same, and the total thickness of heat pellets and upper and lower insulation layers are constant, that is: where d at is the total thickness of the heat pellets and the insulation layers on the upper and lower sides, d at = 51.  From Jeong et al. [14], it is known that the temperature of the 13th layer electrolyte decreases the fastest. In other words, the discharge time of the thermal battery is determined by the temperature of the 13th electrolyte layer. In order to extend the discharge time of the thermal battery, we can increase the thickness of the bottom insulation layer or the thickness of the bottom heat pellet. However, the battery temperature may be too high when the thickness of the bottom heat pellet is increased, resulting in the battery failing. Therefore, we increase dins2 by reducing da and db to achieve longer discharge time of the battery, namely: where Δda and Δdb are the decrease of the thickness of group a and b respectively, and Δdins2 is the thickness difference of the insulation layer at the bottom.
For a thermal battery to reach full discharge phase, the temperature of the electrolyte should always be above its melting point. In order to study the effect of the thickness of the heat pellets on the temperature of the electrolyte, simulations are performed for the thermal battery with different Δda and Δdb. Figure 5 is the temperature profile of the inner part of the thermal battery module at t = 1250 s for different thicknesses, which shows the temperature of the mid-point of the thirteenth layer (the From Jeong et al. [14], it is known that the temperature of the 13th layer electrolyte decreases the fastest. In other words, the discharge time of the thermal battery is determined by the temperature of the 13th electrolyte layer. In order to extend the discharge time of the thermal battery, we can increase the thickness of the bottom insulation layer or the thickness of the bottom heat pellet. However, the battery temperature may be too high when the thickness of the bottom heat pellet is increased, resulting in the battery failing. Therefore, we increase d ins2 by reducing d a and d b to achieve longer discharge time of the battery, namely: where ∆d a and ∆d b are the decrease of the thickness of group a and b respectively, and ∆d ins2 is the thickness difference of the insulation layer at the bottom. For a thermal battery to reach full discharge phase, the temperature of the electrolyte should always be above its melting point. In order to study the effect of the thickness of the heat pellets on the temperature of the electrolyte, simulations are performed for the thermal battery with different ∆d a and ∆d b . Figure 5 is the temperature profile of the inner part of the thermal battery module at t = 1250 s for different thicknesses, which shows the temperature of the mid-point of the thirteenth layer (the layer has the lowest temperature). It can be seen that the internal temperature at the middle of the battery is high while that at both sides is low. From Figure 5a-c, with increasing ∆d b , the temperature of the 13th layer electrolyte increases, which indicates that reducing the thickness of the heat pellets of group b may extend the discharge time of the thermal battery. When ∆d b = 0.24 mm, the temperature of electrolyte increases obviously. This is because the heat pellets become thinner leading to less heat released. However, the heat loss is reduced when the thickness of the insulation is increased, while resulting in a higher battery temperature. It is determined by the heat production of heat pellets and the insulation effect of insulation layers. When ∆d a = 0.36 mm, the mechanism is just the opposite when the heat production of heat pellets plays a decisive role.  From Figure 5a,d, it can be seen that the temperature of the 13th layer electrolyte can be increased by reducing the thickness of the heat pellets of group a. Comparing Figure 5a,e, we find that the temperature of the 13th layer electrolyte is actually decreased when the thickness of the heat pellets in group a is reduced by 0.24 mm. This is because heat pellets become thinner, resulting in less heat generated. At this time, the insulation layer's heat preservation effect is insufficient to compensate for the heat loss caused by the thinning of the heat pellets. Figure 6a shows the influence of the thickness of group a for heat pellets on the time corresponding to the start of solidification of each layer electrolyte. The electrolyte on both sides of the electroactive cells is solidified before the middle part, and the 13th layer electrolyte is solidified firstly, which exhibits that the discharge time of the thermal battery is determined by the solidification time of the 13th layer electrolyte, that is t1, 13 = te-d. When Δda = 0.12 mm, the solidification time of the 13th layer electrolyte is 1290 s, which indicating that the discharge time of the thermal battery is 1290 s. Compared with the results given in Ref. [14] (the discharge time of the thermal battery is 1250 s), the discharge time of this work is increased by 40 s. When Δda = 0.24 mm, the discharge time is 1257 s. It also can be observed that the discharge time is nonlinear with Δda. From Figure 5a,d, it can be seen that the temperature of the 13th layer electrolyte can be increased by reducing the thickness of the heat pellets of group a. Comparing Figure 5a,e, we find that the temperature of the 13th layer electrolyte is actually decreased when the thickness of the heat pellets in group a is reduced by 0.24 mm. This is because heat pellets become thinner, resulting in less heat generated. At this time, the insulation layer's heat preservation effect is insufficient to compensate for the heat loss caused by the thinning of the heat pellets. Figure 6a shows the influence of the thickness of group a for heat pellets on the time corresponding to the start of solidification of each layer electrolyte. The electrolyte on both sides of the electroactive cells is solidified before the middle part, and the 13th layer electrolyte is solidified firstly, which exhibits Energies 2020, 13, 6477 8 of 17 that the discharge time of the thermal battery is determined by the solidification time of the 13th layer electrolyte, that is t 1,13 = t e-d . When ∆d a = 0.12 mm, the solidification time of the 13th layer electrolyte is 1290 s, which indicating that the discharge time of the thermal battery is 1290 s. Compared with the results given in Ref. [14] (the discharge time of the thermal battery is 1250 s), the discharge time of this work is increased by 40 s. When ∆d a = 0.24 mm, the discharge time is 1257 s. It also can be observed that the discharge time is nonlinear with ∆d a . This can be explained that the thickness of the heat pellets becomes thinner leading to a thicker insulation layer which slows the heat exchange between the battery and the air. Therefore, the time for the electrolyte to solidify is prolonged. Figure 6b shows the influence of the thickness of group b heat pellets on the time. It can be seen that when ∆d b = 0.24 mm, the solidification time of the 13th layer electrolyte is 1293 s, namely, the discharge time increases by 43 s. When ∆d b = 0.36 mm, the discharge time is 1265 s, which shows that the discharge time is also nonlinear with ∆d b .
Energies 2020, 13, x FOR PEER REVIEW 9 of 18 pellet. As a result, the heat absorbed by the unit cell decreases, the rate of electrolyte melting decreases, and the activation time increases. Figure 6d shows the influence of the thickness of group b heat pellets on the time. It can be seen that when Δdb = 0.24 mm, the electrolytes of layers four to nine melt almost simultaneously, and the time for complete melting is 1.078 s, that is to say, the activation time is increased by 0.02 s. It can also be seen that the larger the Δdb is, the longer the activation time of the thermal battery is. This is because the heat pellets adjacent to these electrolytes become thinner, resulting in a decrease in heat absorbed by the electrolytes, which causes an increase in activation time.
In the actual design process of thermal battery, the activation time should be as short as possible while the discharge time should be as long as possible. Therefore, it is necessary to comprehensively consider the relationship between the activation time and the discharge time to get the optimized Δda and Δdb to meet the needs of the project. The specific optimization details are given in Section 5.

The Theoretical Basis of Genetic Algorithm
Based on MATLAB function fitting toolbox and genetic algorithm toolbox, we use Δda and Δdb to undertake independent variables. Moreover, the discharge time t1 and activation time t2 are  Figure 6c shows the effect of the thickness of group a heat pellets on the melting time for all layers of electrolyte. It can be seen that the 2nd layer electrolyte finally melts when the time is 1.073 s for ∆d a = 0.12 mm. This can be considered that the activation time of the thermal battery is 1.073 s, namely, t acti = t 2,2 = 1.073 s. The activation time is longer than that given in Ref. [14] (the activation time of the thermal battery is 1.058 s). This indicates that reducing the thickness of the heater can increase the activation time It can also be seen that the larger ∆d a is, the longer activation time of the thermal battery is. This is because the total heat is released when decreases the thickness of the heat pellet. As a result, Energies 2020, 13, 6477 9 of 17 the heat absorbed by the unit cell decreases, the rate of electrolyte melting decreases, and the activation time increases. Figure 6d shows the influence of the thickness of group b heat pellets on the time. It can be seen that when ∆d b = 0.24 mm, the electrolytes of layers four to nine melt almost simultaneously, and the time for complete melting is 1.078 s, that is to say, the activation time is increased by 0.02 s. It can also be seen that the larger the ∆d b is, the longer the activation time of the thermal battery is. This is because the heat pellets adjacent to these electrolytes become thinner, resulting in a decrease in heat absorbed by the electrolytes, which causes an increase in activation time.
In the actual design process of thermal battery, the activation time should be as short as possible while the discharge time should be as long as possible. Therefore, it is necessary to comprehensively consider the relationship between the activation time and the discharge time to get the optimized ∆d a and ∆d b to meet the needs of the project. The specific optimization details are given in Section 5.

The Theoretical Basis of Genetic Algorithm
Based on MATLAB function fitting toolbox and genetic algorithm toolbox, we use ∆d a and ∆d b to undertake independent variables. Moreover, the discharge time t 1 and activation time t 2 are chosen as the objective functions in the process of the genetic algorithm. The flowcharts of algorithm of fitting function and GA technique are shown in Figure 7.

The Generation of the Fitting Function
Two design variables are used to derive the expression of tacti and te-d for Δda = 0, 0.12, 0.16, 0.2, 0.24 and 0.32 mm when Δdb equals 0, 0.18, 0.24, 0.3 and 0.36 mm, respectively. The most common method to solve the curve fitting problem is the least square method. By observing the simulation data (shown in Figure 8) and referring to the simulation work of ref. [18], we find that the change law of these data satisfies the polynomial fitting type. In addition, the R-square (R 2 = 0.9681 and 0.9929) and root mean square error (RMES = 0.0074 and 0.0011) of the fitting function meet the requirements of high-quality fitting, which further proves that the selected fitting type is appropriate. Rewritten Δda/H as x1, Δdb/H as x2 when H = 1 mm, we have the following expressions:

The Generation of the Fitting Function
Two design variables are used to derive the expression of t acti and t e-d for ∆d a = 0, 0.12, 0.16, 0.2, 0.24 and 0.32 mm when ∆d b equals 0, 0.18, 0.24, 0.3 and 0.36 mm, respectively. The most common method to solve the curve fitting problem is the least square method. By observing the simulation data (shown in Figure 8) and referring to the simulation work of ref. [18], we find that the change law of these data satisfies the polynomial fitting type. In addition, the R-square (R 2 = 0.9681 and 0.9929) and root mean square error (RMES = 0.0074 and 0.0011) of the fitting function meet the requirements of high-quality fitting, which further proves that the selected fitting type is appropriate. Rewritten ∆d a /H as x 1 , ∆d b /H as x 2 when H = 1 mm, we have the following expressions: where t e-d, x 1 =x 2 =0 is the discharge time of the thermal battery when x 1 = x 2 = 0, while t acti, x 1 =x 2 =0 is the activation time when x 1 = x 2 = 0. The value of t e-d, x 1 =x 2 =0 is 1255s and that of t acti, where Ni is the control volume output, i N is the average of control volume output, Pi is the predicted value, and n is the total number of data. The comparison of the predicted y1 and y2 with numerical results are shown in Figure 8a,b. The value of R 2 for y1 is 0.9681, and that of RMSE and U2 is 0.0074 and 0.00756 respectively. R 2 = 0.9929, RMSE = 0.0011 and U2 = 0.00107 for y2, which reveals good agreements between the predicted values and the CFD data.

Multi-Objective Optimization
From the above study, it can be found that reducing the thickness of the heat pellets can prolong the discharge time of the thermal battery, which inevitably leads to an increase in the activation time. Therefore, multi-objective optimization of the activation time and discharge time is needed. According to Equations (11) and (12), we need find the maximum value of y1 and the minimum value of y2. Therefore, consider the following objective function: R-square (adj. R 2 ), root mean square error (RMSE) and Theil-U factor (U 2 ) [22] are used to determine the accuracy of the above polynomials and to predict the discharge time and activation time of the thermal battery: where N i is the control volume output, N i is the average of control volume output, P i is the predicted value, and n is the total number of data. The comparison of the predicted y 1 and y 2 with numerical results are shown in Figure 8a,b. The value of R 2 for y 1 is 0.9681, and that of RMSE and U 2 is 0.0074 and 0.00756 respectively. R 2 = 0.9929, RMSE = 0.0011 and U 2 = 0.00107 for y 2 , which reveals good agreements between the predicted values and the CFD data.

Multi-Objective Optimization
From the above study, it can be found that reducing the thickness of the heat pellets can prolong the discharge time of the thermal battery, which inevitably leads to an increase in the activation time. Therefore, multi-objective optimization of the activation time and discharge time is needed. According to Equations (11) and (12), we need find the maximum value of y 1 and the minimum value of y 2 . Therefore, consider the following objective function: It is very important for industrial production to choose the best solution to meet the engineering requirements (i.e., shorten the activation time and extend the discharge time) from the optimization results.
As solutions shown in Pareto fronts are non-dominated, it is necessary to selected the best solutions from Pareto fronts according to their actual applications. TOPSIS [23][24][25] method is adopted to select the optimal solution. The basic principle is that if the solution is the closest to the positive ideal solution (i.e., solution with the minimum 1/y 1 and y 2 ) and the furthest away from the negative ideal solution (i.e., solution with the maximum 1/y 1 and y 2 ), it is the best. Otherwise, it is not optimal. The primary steps are as follows [26]: (1) Create a normalized decision-making matrix (t ij ) m×n by vector normalization: where m and n are the number of solutions and objective functions, respectively. x ij is the element of matrix. (2) Constitute a weighted normalized matrix: where w j is the weighting coefficient. (3) Determine the positive and negative ideal solutions, A + and A − : A + = (min[a 11 , . . . , a m1 ], min[a 12 , . . . , a m2 ], . . . ,min[a 1n , . . . , a mn ]), (19) A − = (max[a 11 , . . . , a m1 ], max[a 12 , . . . , a m2 ], . . . ,max[a 1n , . . . , a mn ]). (20) (4) Calculate the distance of each non-dominated solution to the positive and negative ideal solutions: (5) Calculate the closeness c i of each non-dominated solution to the positive ideal solution: (6) Rank the c i , and the best solution is In the MATLAB toolbox, we use the gamultiobj solver to optimize the expression. In the multi-objective optimization algorithm, the size of population (P) and generation (G) are key parameters affecting the obtained Pareto front. Sensitivity analysis of these two parameters is carried out. Figure 9 is the Pareto front obtained by gamultiobj with different population and generations. It can be seen that with the increase of generation, the Pareto front gradually converges to the same position. The greater the population, the more extensive the distribution of Pareto front. Considering the influences of the value of population and generations, the trial with parameters P = 50, G = 400 has the best performance among all trials, which is used as the obtained Pareto front as showed in Figure 10. The parameters are shown in Table 2. In order to verify the stability of the optimization program, the program is run several times. Figure 10 shows the multi-objective optimization results of four runs. The running results are almost consistent, which proves the stability of the optimized program. It can be seen that y 2 is negatively correlated with the reciprocal of y 1 . That is to say, the longer activation time is, the longer discharge time is. However, we need reduce the activation time and increase the discharge time. The points in the graph are all non-inferior solutions, which has smaller target conflicts than other solutions. Table 3 shows all non-inferior solutions and structural parameters for the fourth run. When the activation time is not more than 1.075 s, we select x 1 = 0 and x 2 = 0.161, that is, ∆d a = 0, ∆d b = 0.161 mm, the optimization result is the best.