Optimal Design of Transmission Shafts: a Continuous Genetic Algorithm Approach

This paper presents an analysis of the optimal design of transmission shafts by adopting the approach of a novel continuous genetic algorithm. The optimization case study is formulated as a single-objective optimization problem whose objective function is the minimization of the total weight that results from the sum of all the sections in the shaft. Additionally, mechanical stresses and constructive characteristics are considered constraints in this case. The proposed optimization model corresponds to a nonlinear non-convex optimization problem which is numerically solved with a continuous variant of genetic algorithms. SKYCIV®and Autodesk Inventor®were used to verify the quality and robustness of the numerical results in this paper by means of simulation tools and analysis. The results obtained demonstrates that the methodology proposed reduce the complexity and improving the results obtained in comparison to conventional mechanical design.


Introduction
The drive shafts of industrial machinery are responsible for transmitting power, generally from the engine to the mechanism, which makes the equipment operate [1]. Several methods are used to improve them. Most of them are iterative and entail long hours of work to develop a single piece because it should meet both physical and minimum safety criteria. This process starts with mathematical modeling, which is solved by trial and error to obtain a design close to the expectations [2]. Other tools that support this process are prototyping and simulation, which validate data and certify that the element under development will withstand the loads it will experience [3]. Different studies have been conducted in this field. For instance, in 2016, Reddy and Nagaraju studied the minimization of the weight of a drive shaft of an automobile [4]. To achieve their goal, they created a mathematical model considering the constructive and physical parameters that must be respected for a correct operation. Subsequently, they carried out a validation with Finite Element Analysis (FEA) software to guarantee compliance with the characteristics mentioned above. Other works of this type [5] have experimentally and computationally evaluated the strength of the joint between a pinion and an engine to calculate different phenomena that affect the coupling. Due to the time and economic resources required by traditional design techniques, optimization algorithms have been introduced. They explore a search space under given parameters and a set of constraints in order to provide an optimal solution to the mathematical model of the problem [6] in a quick and accurate manner. Thus, they add the value of safety because, if the mathematical model is analyzed with specialized computer programs, the human 803 error that occurs when the conventional method is applied would be reduced [7]. Several algorithms have been implemented in mechanical design to solve problems such as the minimization of the weight of a cantilever beam, the minimization of the cost of storage tanks, the reduction of the weight of supports, and the minimization of the cost of speed bumps for roads, among others. These problems have been solved by means of the Particle Swarm Optimization (PSO) algorithm, League Championship Algorithm (LCA), and Multi-Objective Differential Evolutionary (MODE) algorithms [7], [9] and [10]. These algorithms contain a set of basic constraints defined by physical design criteria and aim at obtaining reasonable and feasible solutions. Genetic algorithms (GA) have also been used for that purpose. For example, by means of a GA, some authors studied machining optimization through Computer Numerical Control (CNC) and geometry variations [11]. Another study [12] used a GA technique to optimize wind turbines for wind farms in order to find the rotor that would produce most power. Since optimization techniques, especially metaheuristic ones [13], are widely used to find a solution to different nonlinear problems [14], this work proposes the minimization of the weight of a drive shaft with an abrupt change in the cross section. In this case, the objective function is the equation to calculate the weight of a shaft, and the variables proposed for the problem are the diameters assigned to each of its sections. The constraints of the problem are the general equation to design a fatigue-safe drive shaft [15], a constructive measure that limits the minimum distance of the steps for the coupling of the mechanical elements and fasteners (established by the manufacturers of bearings [16]) and the maximum allowable deflection, which is defined by calculating the elastic curve [17]. In order to solve the case described above, this study proposes a GA because it is a common technique to address this type of problems and it provides excellent results in multi-variable optimization with a low computational cost [18]. This document is divided as follows. Section 2 presents the analysis of the problem, type of shaft, and the elements coupled to the latter. Section 3 describes the mathematical model proposed to calculate the dimensions of the shaft. Said model is composed of a single-objective function and the set of restrictions that represent the problem. Section 4 explains the GA proposed to solve the mathematical model described in the previous section and also presents the method proposed to solve the mathematical model based on the parameters required to implement the algorithm; besides, it details the way the load analysis and the validation were carried out. The results obtained by the GA are found in Section 5, along with the corresponding discussion and verification by means of open-source software for beam calculation by SKYCIV®and Autodesk Inventor®. Finally, Section 6 includes the conclusions of this study.

Problem Formulation
Essential elements in any industrial machine that performs a mechanical movement, drive shafts may be subject to loads that change over time. As a result, they should exhibit physical characteristics of material strength and fundamental design criteria that enable their correct operation without presenting premature failures. Due to its importance in the industry [15], this work examines a stepped drive shaft that experiences a concentration of stresses caused by abrupt changes in the cross section of each step. The main objective of the proposed design approach is the minimization of the weight of the shaft while respecting a set of constraints that represent that part. The study case in this article was taken from another work [19] that describes a stepped drive shaft divided into 3 transversal sections, whose lengths are L1 = 5.9in, L2 = 19in, and L3 = 1in.
Furthermore, the physical analysis considers a group of different components and their characteristics. A pulley of 20in pitch diameter and a weight of 60lbs. is coupled to point A; said pulley reduces the shaft speed. Besides, a spur gear of 10in in pitch diameter and 25lbs. is coupled to point C. Pillow block bearings are coupled to support points B and D. The transmitted power is 20hp. The drive shaft is manufactured with AISI 1040 CD steel, and the assigned factor of safety is 2.2. Figure 1 presents an illustration of the location of the points and the elements, where the former have the following values: A = 0in, B = 5.9in, C = 15.9in, and D = 25.9in.
In order to calculate the optimal diameters of the shaft, this work proposes the use of a GA, a meta-heuristic optimization method [20]. The GA will explore the solution space of the problem while evaluating, at each iteration, the weight function and the constraints described in the mathematical model of the problem.

Mathematical Modeling
The main objective of the mathematical model of the problem described in this study is to define and guarantee that the diameters selected by the algorithm will be adequate and, at the same time, that they will represent a reasonable and feasible solution to manufacture the shaft. For that reason, the algorithm must include formulas that describe the problem and conduct analyses based on statics, dynamics, the mechanics of flexible materials, and mechanical design methods.
Equation (1) represents the adaptation function that corresponds to the reduction of the weight of the shaft as a function diameters d 1 , d 2 and d 3 , assigned to each segment of the shaft, L 1 , L 2 and L 3 , respectively.
where γ is the specific weight of AISI 1040 CD steel, 0.2834lb/in 3 . This problem is limited by a set of constraints that were considered in the formulation of the mathematical model of the problem. They result from the physical analysis of the element under study, which was conducted during the design stage. The main constraint corresponds to the general equation to design a fatigue-safe drive shaft [15] as a function of the diameters selected by the designer or the solution algorithm. Such general equation defines the minimum diameter the shaft should have to withstand the loads it will experience. Equation (2) presents this constraint in a general form.
where N f is the desired factor of security, 2.2; S ut , the maximum resistance to mechanical stress, 75000 psi; and S f , fatigue strength corrected in a selected life cycle. The latter is calculated by applying Norton's theory, which considers corrections to ultimate strength experimentally reported for this material in websites such as MatWeb [21]. The value of fatigue strength in this case was 24314.3354 d psi. K f is the stress concentration factor for the alternating component of the bending stress (2.05, 3, and 2.4, for d 1 , d 2 and d 3 , respectively). K f s denotes the stress concentration factor for the alternating component of the torsional shear stress (2.05 for the three diameters). K f m represents the stress concentration factor for the mean normal component (which has the same values as K f ). Finally, K f sm is the stress concentration factor for the shear mean component (2.05). These concentration factors were taken from the theories by [15] and [1]. M a and T a refer to alternating moment and torque; M m and T m , to average moment and torque, which are calculated by analyzing bending loads and applying statics theory [22]. Subsequently, a free body diagram of the shaft in each orthogonal plane was drawn to calculate maximum and minimum moments. Besides, a 2D vector analysis was conducted in the "horizontal" and "vertical" planes considering the maximum and minimum forces. The tool SKYCIV® [23] for beam calculation was used to draw internal shear and bending stress diagrams that served as a basis to calculate, by vector addition, the critical loads the drive shaft experienced. Subsequently, the state of the load was characterized by calculating average and alternating moments. Since every segment of the drive shaft undergoes different forces, the alternating moments that occur in B and C were calculated (1154.7 and 1097.7 lbin, respectively) as well as their average moments (1619.8 and 1828.8 lbin, respectively). Furthermore, the values of average and alternating torque of d 1 and d 2 , reached 2187.5 and 1312.5 lbin, respectively. Torque is not considered for d 3 because this section does not experience that force.
The constraint of meeting the minimum diameters is expressed in equation (3), where d i represents the diameter of each cross section under analysis at each iteration. Said equation establishes the difference between the selected diameter and the minimum diameter analyzed in equation (2), whose difference should be equal to or greater than 0. Since this drive shaft is composed of 3 cross sections, a constraint should be defined for the diameter applied to each section. The resulting constraints are presented in equations (4) to (6).
Additionally, a restriction that determines the step between the sections of the shaft should be defined: minimum 0.0787 in with respect to the diameter of section L 2 . This limitation is derived from constructive criteria established by bearing manufacturers [16] (equation (7) and (8)).
Equation (9) includes a constraint associated with the maximum allowable deflection in the shaft, which cannot be greater than 0.005 in. This step is taken to guarantee that the gears operate properly because, if such value of lateral movement is exceeded, the elements could come apart and thus stop transmitting or transmit incorrectly, which would later cause system failures [15].

Proposed Method
To solve the mathematical model of the problem, this work proposes a GA, which is explained in the next section. Furthermore, the processes of load analysis and verification by software are also described.

Genetic Algorithm
GA is a classical optimization technique [24], have been used multiple times to solve continuous optimization problems (e.g. optimization of nonlinear continuous functions and second-order borderline differential equations) with satisfactory results. The following are the 5 characteristics of a GA: (a) Generation of the initial population. Each one of these characteristics is greatly important to satisfactorily solve any optimization problem with said tool. For that reason, each one of them is expanded in the following section.

Initial Population
This is the first step in any optimization technique. For this problem, this study proposes a population with a size of a rows and s columns (i.e. an axs matrix) that uses the code in Figure 2. Where a corresponds to the number of potential solutions (called group of individuals) and s is the number of variables in the problem. In this case, since the shaft has three sections, 3 columns are considered. This initial population will be calculated as shown in Table 1. Where d represents each value of diameter assigned by the algorithm to each position of the matrix of individuals. This value is initially calculated at random, respecting the minimum and maximum limits established for the diameters.
Moreover, the way the problem will be coded should be described in order to define the individual values that will be entered into the model analyzed at each iteration. A 1-row 3-column (1x3) vector was used for this problem; its columns contain the three diameters assigned to lengths 1, 2, and 3, in that order. The values this vector takes depend on the set of constraints presented in Section 4 ( Figure (2)). Section

Objective Function
In meta-heuristics theory, the objective function corresponds to the performance of the function assigned to any individual contained in the population, i.e. said function determines the quality of the arbitrary solution. It should be noted that GAs solve optimization problems by transforming a restricted problem to a conditional one. For that reason, the objective function of this problem will be weight minimization in addition to meeting the set of constraints that will take values depending on the breach of certain conditions under analysis. This part of the process aims at being able to classify a configuration of diameters that do not meet one of the conditions restricting the problem as a good solution. This is expressed in equation (10).
P en = (p 1 + p 2 + p 3 + p 4 + p 5 + p 6 )F p Equation (11) represents the penalties associated with the set of constraints that represent the problem, which is necessary to enable the optimization technique to advance on the infeasibility region and, at the same time, improve the exploration of the search space. Said equation is expressed as the sum of penalties multiplied by a penalty factor (F p = 1.5e3, found by trial and error), which enables the sum of Penn with f1, thus penalizing the solutions that violate any constraint. Finally, the identification and interpretation of each one of the constraints that define the problem resulted in the set of penalties employed in equation (11), which is composed of equations (12) to (17). Such penalties are associated with each one of the restrictions previously mentioned and defined by the maximum value function, where, after each constraint is analyzed, the function will take the maximum value as a result, i.e. a value of g if there is violation or 0 if it meets the criteria.

Offspring
Since a GA is an iterative process, new potential solutions to the problem under study should be generated in order to replace bad solutions contained in the current population. To create this new group of individuals, the operators of mutation and recombination are adapted. In addition, classical selection is used to enable good solutions to stay in the group of individuals, thus allowing them to improve their position [25].
i Selection: The offspring arbitrarily select a subgroup of individuals contained in the current population. During this step, a random number r from 1 to a is chosen, i.e. r = (a -1) rand. If r < a, an additional (ar) x s matrix is generated with new potential solutions using the same strategy as with the initial population. The total group of selected individuals is defined by the two strategies mentioned above. ii Recombination: This process alters the offspring that follows the initial population. If the probability of recombination r p is higher than 50%, two arbitrary individuals (randomly selected) are averaged to generate a new potential individual. It should be noted that this operation always produces feasible individuals from the initial population, and they are as good as the random solutions generated within the admissibility region of the limits of the problem. This process will continue to obtain an offspring with a potential solutions. iii Mutation: This step explores the probability of m p mutation. In other words, if m p is over 50% (a value that has been arbitrarily selected), a new random position of the potential solution is modified by an arbitrary diameter value that ensures compliance. If m p is under 50%, the potential solution is not modified. This process is applied to all the individuals in the offspring.
Once the offspring has been generated, the objective function is evaluated again.

New Population
The new population will save the group of best solutions found by the GA until the current iteration. Afterward, a new population is generated by the combination of the group of individuals in the offspring, which produces a population with 2a potential solutions. As a result, two potential solutions are identical; hence, one of them is eliminated from the list. This process is repeated until it can be guaranteed that all the potential solutions are different. Subsequently, the resulting list of potential solutions is organized in ascending order, with all the individuals as a function of their objective function. Besides, the first a potential solutions are selected as a new population that will go on to the next iterative cycle.

Stopping Criterion
The continuous GA proposed in this work stops the optimization process when some of the following conditions are satisfied: i The total number of iterations was reached.
ii The best potential solution does not improve after m continuous iterative cycles.
Otherwise, the algorithm will continue generating offspring nonstop. The algorithm 1 shows the pseudo-code that represents the GA employed to solve the problem under analysis.
Data: genetic algorithm parameters, parameters of the mathematical model.
Result: Impress results Break; end end Algorithm 1: Proposed pseudo-code for the Genetic-Algorithm for solving the optimal diameters of the drive shaft.

Parameters of the GA
To be implemented, the algorithm requires a series of parameters, which are listed in Table 2. These values should be entered into the algorithm, since they will enable it to correctly advance over the solution space.

Load Analysis
This analysis considered the torque, calculated based on the angular velocity and the power of the system in each scenario (maximum and minimum load). Afterward, average and alternating torque were computed, which enabled to estimate the forces exerted by the pinion and the pulley on the shaft in order to subsequently decompose said forces into their vertical and horizontal components. Online software for beam calculation, SKYCIV®, was used to verify the internal bending moment the shaft undergoes due to external forces that act upon the places where transmission elements and bearings are located. Figures 3, 4, 5 and 6 present the free body and bending moment diagrams of the vertical and horizontal planes of the shaft in both scenarios. Based on this information, maximum and minimum moments are identified. They will enable to calculate the alternating and average moments used in equation (5). Additionally, it can be observed that C is the critical point in both scenarios.

Deflection Analysis
This analysis was carried out using a mathematical model that describes the elastic curve [17] (Equation 18).
Where F a , F b , and F c correspond to the value of the forces in the vertical plane and the value of the singularity functions changes depending on the point to be analyzed. Such value becomes 0 if the result of the subtraction between the distance under analysis and the location of the force is a negative number. The constants are the result of the double integration of the singularity equation that describes the moment (19). E refers to Young's modulus (29000 ksi for this material) and I represents the geometric property of the area.

Verification with SKYCIV®
The software for beam calculation provided by SKYCIV® enables to calculate the elastic curve based on the configuration of the length of the element, the loads that act upon it, the supports in place, the elastic modulus, and

Verification with Autodesk Inventor®
The shaft generator provided by Autodesk Inventor® was employed to conduct a more thorough verification of the deflection. Such tool enables to carry out a FEA simulation of a shaft and to enter values of diameters, lengths, and radii in accordance with the values and the location of the forces that act upon it, the location of the bearings, and the properties of the material under analysis. The next step is to calculate the maximum deflection the shaft experiences along its span. Figure 8 shows the module of the software that generates the shafts.

Results and Discussion
This work implemented an optimization GA to solve the model that enables to calculate the diameters of the cross sections of a shaft under the conditions studied in this article. The following section presents the results obtained by the algorithm programmed in Matlab® on a HP Z600 desktop with 12Gb of RAM memory and 6 cores. Such solution was validated in SKYCIV® and Autodesk Inventor®, which are used to verify both the elastic curve and maximum deflection.

Results of the GA
The computation time it takes the genetic algorithm to solve the model proposed in this study is less than 3s. Table 3 presents the best configuration of diameters found by said algorithm. After being analyzed, this solution is considered valid because it does not violate any of the constraints described in the Method section. This proves that a drive shaft with the calculated diameters will withstand the loads it will experience. Besides, it can be noticed that the difference between d 1 and d 2 with respect to d 3 is higher than 0.0787 in and lower than 0.5 in, which represents an adequate step to couple the bearings. Additionally, the deflection analysis estimated that the shaft would exhibit a maximum displacement of 0.001693 in, which is under the maximum value, 0.005 in. Therefore, the parts would satisfy the design criteria for rigid drive shafts. With this configuration of diameters, the shaft would weigh 17.6834 lb, which is the minimum possible weight of the shaft under study that respects the constraints defined in the optimization algorithm. A comparison of these results with the literature reveals a noticeable difference with a particular study [19]. In that case, the minimum weight was 1.688e6 lb, which constitutes an in-feasible solution, physically and construction-wise. After this value was identified and a thorough analysis was conducted, erroneous calculations were found in the constraints as well as the objective function of said work. On one hand, the value of that answer is the result of the fact that the specific weight the authors used in that case does not correspond to the material in the problem description. On the other hand, the set of constraints they propose establishes a variation in the lengths in order to find the optimal diameter. This causes the moments that occur along the shaft to change; thus, they would need to be calculated at each iteration, which would require a more in-depth analysis to produce a correct solution. This point was not considered by the authors of that work.

Verification with Software
The following are the results of the verification with software; they guarantee compliance with the deflection constraint.

SKYCIV®
The open-source software for beam calculation [23] enables to compute moments and evaluate elastic curves. The latter express the maximum deflection depending on the loads experienced by the shaft under analysis, the bearings in place, the moment of inertia, and the Young's modulus of the material. When the value of moment of inertia of the calculated solution (0.5034 in 4 ) is entered, the software generates an elastic curve of the vertical plane of the shaft in conditions of maximum load (Figure 9), which shows the separation between gears during transmission. This elastic curve allows to conclude that, although the greatest deflection is found in point A, it does not represent a cause for penalty because the flexible transmission system composed of pulleys and a belt is not affected by said value. Furthermore, it can be observed that the maximum deflection between points B and D is not located at 15.9 in (point C) but at 12.163 in; however, it does not represent a considerable value. As previously mentioned, the maximum deflection found on point C of the solution is 0.001693 in, which would not affect the performance of the coupled gears.

Autodesk Inventor®
Afterward, specific data was entered into the shaft generator module of Autodesk Inventor®: the diameters found by the GA and the maximum forces the shaft experiences in the vertical plane (where the maximum deflection of point C should not exceed 0.005 in). Figure 10 shows the results of the module with the location of the forces and the bearings. It can be observed that the maximum deflection on point C is 0.001788 in. Besides, the correspondence between the weight of the shaft and the value calculated in the mathematical model was verified. When the analysis was conducted in Matlab®and the mathematical model was solved with the selected diameters, the deflection equation described in the Method section yielded a value of 0.001693 in, which represents a 6% error and the total satisfaction of the design criteria mentioned above.

Conclusions
In this paper a new approach for the continuous genetic algorithm is proposed for mechanical design of transmission shaft. Using a single-objective formulation that to employ the reduction of the weight of the shaft as objective function, and the mechanical technical criteria and the deflection associated to this type of machine element as set of constrains. The combination between the mathematical model proposed and the GA, allow to reduce the complexity and improving the results obtained in comparison to conventional mechanical design.
The GA employed in this study proved to be useful and quick because it found the best configuration of diameters at a low computation cost, which makes it an efficient means to solve the problem discussed in this article.
This strategy avoids conventional iterative work and allows to solve a great number of mechanical design problems. It only requires an adequate mathematical model that describes the physics of the system. Introducing this type of optimization algorithms provides great help to improve or design different elements of industrial machinery while guaranteeing that the solutions are reliable, adequate, and feasible, physically and constructionwise. The verification tool also added value in terms of safety as it demonstrated that the mathematical model was correctly formulated. The error found between the solutions provided by Matlab® and Autodesk Inventor® is acceptable and, since the software uses FEA, it offers mathematical analysis features that enable to produce results closer to real values.