Macro and Micro-Geometric Robust Design of Cylindrical Roller Bearings Using Fatigue and Wear Considerations


 Cylindrical roller bearings have been chosen for the optimization using a robust design because tolerances are always provided during manufacturing. For rolling element bearings, a measure of satisfactory performance is a long life, lubrication and thermal characteristics. Hence, the first two objectives have been optimized simultaneously using the robust design. For bearing design, variables chosen include four basic internal dimensions and five design constraint factors. Genetic algorithm is implemented for solving the problem by considering both the objective functions, simultaneously. The lives of the bearing are found to be improved at least two times that of the standard bearings. It has been found that during manufacturing, designer should have more attention on the roller diameter compared to others dimensions.


Introduction
Rolling bearings are the bearings that is used in rolling elements for carrying the loads. Rolling elements may be the balls or rollers depending on the amount of load the bearing has to carry and the radial and axial space. These are situated between the annular space of two races/rings, which is generally referred as the outer ring and inner ring. There exists a relative motion between these rings that causes the rolling elements to roll with negligible friction. Based on the shapes of the rolling elements, rolling bearings are of two types namely the ball bearings and the roller bearings. Roller bearings have many applications than the ball bearing due to their ability to carry large radial load. And this is due to the fact that the type of contact between the rolling elements and the raceways in roller bearing is line and in ball bearing is point, so a large stressed develops in ball bearing for a same value of load on them. Cylindrical roller bearing (CRBs) is, such an example of roller bearing, which has the above mentioned advantage. In addition to this, CRBs have exceptionally low-friction torque characteristics that make them suitable for high-speed operations. These characteristics enables the CRBs to be used widely in aerospace applications, automobile sectors, railways, wind or gas turbines, etc. Thus, the bearings must be designed such that they will have their maximum performance. This provokes the researchers to develop new technologies that will enhance the bearings performance and hence their lives.
Changsen [1] gave an organized, clear and logical sequence of topics regarding analytical methods of rolling bearings. He systematically explained the geometry, kinematics, load distribution, stress and deformation in rolling bearings for calculation of the static load rating, dynamic load rating, life of rolling bearings and elasto-hydrodynamic lubricant film thickness. His work also gave insight to optimization of design and computer aided design of rolling bearings. Harris [2] analyzed the rolling bearings thoroughly, like their types, structure, tolerance, simple and complex load distribution, stress analysis, thermal analysis, etc. Analytical solutions to various problems have been represented and empirical approaches to certain problems have been used where they have appeared to be more practical. Seireg [3] gave a brief review of the use of optimization methods in the design of various machines parts. Deb et al. [4] developed the non-dominated sorting genetic algorithm (NSGA-II), which is being used in the current work. A fast non-dominated sorting approach with lesser computational complexity than other optimization algorithms is presented. Also, a selection operator is presented that creates a mating pool by combining the parent and offspring populations and selecting the best (with respect to fitness and spread) solutions. Simulation results on different test problems showed that the proposed NSGA-II, in most problems, is able to find much better spread of solutions and better convergence near the true Pareto-optimal front (POF) compared to other elitist multi-objective evolutionary algorithms (MOEAs). Moreover, the definition of dominance has been modified in order to solve constrained multi-objective problems efficiently. Branke et al. [5] discussed a modified MOEA, which is able to focus search on knee regions, resulting in a smaller set of solutions that are likely to be more relevant to the decision maker. Maday [6] used bounded variable methods of the calculus of variable to determine the optimum configuration for one-dimensional hydrodynamic gas slider bearings. The design criterion was maximizing the load carrying capacity of the bearing. Jamison [7] discussed the effects of tolerance control in manufacturing and wear at rib-roller contact under poor lubricating conditions on bearing performance.
Chakraborty et al. [8] described design optimization problem of rolling element bearings with five design parameters using genetic algorithms (GA) based on the requirements of longest fatigue life. The main limitation of the method was that the use of unrealistic constraints. The assembling angles were assumed and values of constraint constants were chosen fixed to solve the optimization problem. Chaturbhuj et al. [9] described design optimization method for tapered roller bearings using GA. Hirani [10] minimized the temperature rise, oil feed flow and power loss in hydrodynamic bearings. Link et al. [11] optimized the suction muffler and bearing by increasing the energy efficiency for the suction muffler and reducing the power losses for bearing using GA. Saruhan et al. [12] used the GA for developing an optimum design method for tilting pad bearings. The approach optimizes the minimum film thickness, power loss and maximum film temperature. Hirani and Suh [13] optimized the fluidfilm static loaded journal bearings by minimizing the oil flow and power loss using a finite difference mass conserving algorithm, Pareto optimal concept, GA and fitness sharing. Wang [14] investigated the performance of parallel optimization by means of the GA for the lubrication analysis. Kim et al. [15] proposed a bearing parameter identification methodology based on global optimization scheme using measured unbalance response of rotor-bearing system. A hybrid evolutionary algorithm, clustering-based hybrid evolutionary algorithm (CHEA), is proposed for global optimization scheme to improve the convergence speed and global search ability. Gupta et al. [16] proposed a design optimization method of rolling bearings with three objective functions, including the static capacity, dynamic capacity and the elasto-hydrodynamic minimum film thickness, using NSGA-II. A parametric sensitivity analysis of the design variables has also been carried out.
Rao et al. [17] described a constraint non-linear optimization procedure based on GAs for designing rolling element bearings based on the maximum fatigue life as objective function. A convergence study has also been performed. Kumar et al. [18,19] used the real coded GA for the optimization of rolling bearings with the objective as maximization of fatigue life. Peng et al. [20] performed numerical calculation on elasto-hydrodynamic lubrication (EHL) of grease lubricated ball bearing, it predicted that the inner ring surface compared with the outer ring is more likely to have wear and fatigue failure. Tiwari et al. [21] proposed an optimization method to design a tapered roller bearing to have an excellent life. A convergence study to check the nature of convergence of design variables and objective function, and a sensitivity analysis to find the variation of objective function with the variation in design variables have also been carried out. Solanki et al. [22] designed the bearings by maximizing their load carrying capacity with respect to clearance, and minimizing the temperature rise and film thickness. Panda et al. [23] considered the contact stress between raceways and rolling elements of ball bearings as a constraint and optimized the bearings using the particle swarm optimization algorithm. Dragoni [24], while optimizing the dimensions of tapered roller bearings for the maximum static capacity, he found that static capacity varies linearly with-roller infill, square of bearing mean/pitch diameter, ratio of roller length to roller diameter. Tiwari et al. [25,26] used the artificial bee colony algorithm for the optimization of rolling bearings. Kim et al. [27] enhanced radial and axial stiffness values of an angular contact ball bearing using a micro-genetic algorithm. Kalyan and Tiwari [28] optimized the design of needle roller bearings for getting an excellent fatigue and wear lives by using NSGA-II. The objective functions were optimized individually and concurrently. Result showed good improvement over standard bearings.
The available literature is focused on the conventional design optimization of rolling bearings with the single and multi-objective functions with GAs (e.g., NSGA-II). However, it is impossible to manufacture bearings with an exact optimal dimension. They may always have some uncertainties, like tolerances while manufacturing, which may be either error or intentionally defined by the designer. No doubt these uncertainties degrade the performance of the bearings. Hence, the design of the bearing should be such that after manufacturing, if it has tolerances then still it will have the same performance as mentioned in the design. This compelled the researchers to think on it. However, Messac and Yahaya in their work [29] partially eradicated this problem. They proposed a flexible physical programming based robust design optimization (RDO) technique, that ensures minimal variability of the objective function in the case of input variations or uncertainties. The method allows the designer to prescribe parameter tolerance levels as well as to specify the desired range of variation in the objective function. Chen et al. [30] explored the effectiveness of the physical programming approach in explicitly addressing the issue of design robustness and established the superiority of physical programming over other conventional methods. It also illustrates that the physical programming method is among the most effective multi-criteria mathematical programming techniques for the generation of Pareto solutions that belong to both convex and non-convex efficient frontiers. Giassi et al. [31] developed a multidisciplinary optimization method based on robust design techniques for application of the concurrent engineering, so as to manage the uncertainty due to design team collaboration. The method assures effective design work distribution. One such major work can also be found in Taguchi's design [32]. The Taguchi's model discussed mainly about loss function. Although Taguchi's method was effective at improving quality, there were several inadequacies in the method, like for the problems showing highly nonlinear the results were unsatisfactory. Verma and Tiwari [33] proposed an optimization method for robust design of tapered roller bearings. Bearings were designed such that they have the maximum possible dynamic capacity along with minimum variation in it, which was due to the tolerances associated with their basic dimensions. They only considered dynamic capacity as the objective function.
The dynamic capacity is not only a factor that affects bearings life and performance, the elastohydrodynamic minimum film thickness is also equally important. Therefore, it must be optimized in addition to dynamic capacity in order to get a bearing with an excellent life. So, with these objectives in mind, in this study, the robust design optimization of CRBs is done using the NSGA-II algorithm. Hence, for the present study, the robust design has been achieved using flexible programming approach, which treats each objective function and its desired tolerance, and the design variables and their maximum possible tolerances independently, while expressing designer preference and hence, eliminating all the above stated deficiencies.

Geometrical aspects of CRBs
A CRB is designated by the standard dimensions-bore/inner diameter (d), outer diameter (D), and bearing width (B). The internal dimensions that specify a CRB completely are pitch diameter (Dm), roller diameter (Dr), roller effective length (le) and number of rollers (Z). The cross-section view of a CRB can be depicted from Fig. 1, which shows various geometrical parameters along with the raceways and guide ring.

Problem construction for robust optimization of CRBs design
Any optimization problem has generally four components namely the design variables/factors, objective functions, constraints and design variable bounds.

Objective functions:
The aim of the designer, in the form of real valued function, which is either to be minimized or maximized is known as objective function. A very known principle called the duality principle that converts the minimization of fi(x) to maximization of "fi(x)" and vice-versa.

Minimize/Maximize
where x refers to the design variable vector that are to be optimized, where, KDmin and KDmax are the parameters that limits the minimum and maximum roller diameter, respectively, e guarantees the kinetic motion between the roller and raceway, ε restricts the thickness of outer ring and β restricts the roller effective length. While the parameters Dr, Dm, le and Z are the internal dimensions of the bearing, the remaining five parameters present in the constraints and provide a feasible design space.
The objective function for the robust optimization will be of the form where ∆f is the deviation in the objective function and can be expressed as It must be noted that  may be the variation in design variable, constraint or objective function. The vector ∆x comprises of where ∆Dm is the pitch diameter tolerance, ∆Dr is the roller diameter tolerance, ∆le is the roller length tolerance. The value of each of these tolerances has to be input by the user.
The present paper considers two objective functions-dynamic capacity (Cd) and elasto-hydrodynamic minimum film thickness (hmin.) that are supposed to be optimized simultaneously.
Dynamic capacity (Cd): Herein, the dynamic capacity is optimized, which is defined as "constant radial load, which can be endured by a group of evidently alike bearings for a rating life of one million revolutions of the inner ring for a static load and fixed outer ring". It can be expressed for a CRB as [1]   29 As we already know, that the objective function varies only due to the design variables that do have tolerances and are present in the mathematical formula of the objective function, hence, in the case of dynamic capacity, the variation will not be due to the bearing standard dimensions D, d, and B. Hence, using Eq. (4), we have Applying log to Eq. (13) and then differentiating with respect to Dr, we get On differentiation of Eq. (12) with respect to Dr and dividing both sides by K, we get Hence, from Eqs. (14) and (16) Similarly, applying differentiation to Eq. (6) with respect to Dm and le, and continuing with the same approach as Eqs.

Elasto-hydrodynamic minimum film thickness (hmin ):
Other than the fatigue life of bearing, the bearing must have good wear life. This wear life depends on the elastohydrodynamic lubrication (EHL) [1]. Across the raceways (both inner and outer), there exist a minimum film thickness (hmin), that decides the bearing's wear life as it avoids the direct metal-to-metal contact between the raceway and the rolling element. As the bearing demands a longer wear life, therefore this hmin should be maximized. Since, there are two raceways (inner and outer) in rolling bearings, therefore individual hmin will exist across the raceways, namely hmin, inner for the inner raceway and hmin, outer for the outer raceway. But this is not necessary that each hmin need to be maximize separately. Only minimum of these two need to be maximize. The mathematical formula for hmin available in [1] is valid for both inner and outer raceways that can be given here as ; For the deterministic design optimization, the aim is to maximize the lesser one among hinner and houter, and can be mathematically written as Using Eqs. (20) through (24) for the inner raceway, we get   Likewise, differentiating Eq. (29) with respect to Dm and le, we get where xl represents the l th element of the vector. However, for the robust optimization, it is expressed as where ∆xl represents the tolerance in the design variable corresponding to the l th element of the vector. It must have noted that tolerance is added to the lower bound and subtracted to the upper bound. This is because the new range that we obtained is still lying within the range stated in deterministic optimization. Hence ensuring the design feasibility along with its robustness.

Constraints:
Constraints help in defining a feasible domain of values for design parameters and variation in the objective functions. Constraints for the deterministic optimization is expressed as while, x still holds the same meaning as explained above, p refers to vector of those parameters in the constraints, that are not present in the objective function, and whose value remain constant during optimization (but may undergo variations during manufacturing and hence their variation must be considered, which we shall deal with, in the robust optimization), g refers to the design constraints for the deterministic optimization and J is the total number of design constraints. Hence, in case of CRBs, we have The vector ∆p comprises of The changes made in constraints of the deterministic optimization to convert it to a robust optimization, as seen in Eq. (43) are as follows Here, the value for each of these tolerances must be input by the user. The constraints from Eq. (43) can further be stated as (47) where h is the design constraint in inequality form. It must be noted that the number of deterministic optimization constraints is different from the number of robust design constraints. Hence, we must represent them by different values "J" and "K" respectively. In Eqs. (43) and (46), we have shown how to convert each deterministic optimization constraint to the robust optimization constraint. When Eq. (43) is being written, we are still utilizing the deterministic constraint only, hence we have written "j" for this equation. While, for Eq. (46) it is already converted to robust design constraint, hence it is shown by "k". But, it must be noted that for all those constraints that are getting converted from deterministic to robust, the values of "J" and "K" are same. Although the range of values shall differ for "j" and "k" due to the presence of two constraints in the robust optimization, prescribing the desired range of variation in the objective function, which we will see in further sections.
Constraints 1&2: From the geometry of the bearing, it is clear that the bearing pitch diameter should be smaller than the outer diameter and should be greater than the inner diameter. But this does not provide an enough constraint. We must also consider chamfering in the corners of rings to provide them a sufficient space. Thus the constraint is Constraints 3&4: The lower and upper bound of the roller mean diameter (Dr) is decided by the strength and geometric aspects, respectively. The least value of Dr should be such that rollers can easily resist the stress caused by the load applied on the bearing. While the upper limit of Dr is constricted by the boundary dimensions D and d, and the chamfering radii of the corners of rings.
where Qmax is the maximum normal load acting between the roller and the raceway, and max c  is the maximum contact stress on the roller. Qmax can be found out using the Stribeck's expression [2].
where Z is the number of rollers, Fr is the radial load acting on the roller, and α (= 0 0 ) is referred as contact angle, which is the angle between the roller and outer raceway.
Constraints 9&10: The constraints that need to be added for ensuring the running mobility of the rollers and raceways are Constraint 11: To sustain the rollers at the bottom of the outer raceway of the bearing, the thickness of the outer ring should be greater than a particular value equal to Dr [25].
Constraint 12: Another constraint for the outer ring considering the chamfering corners can be added as Constraint 13: From the geometry of the bearing, it can be seen that the inner ring has convex shape, due to which it gets more stress compared to the outer ring. Therefore, to make the inner ring stronger, its thickness must be greater than that of outer ring.
Constraint 14: Due to cylindrical shape of the roller, the type of contact between the roller and the raceways is a line contact. This line contact leads to generation of a very high stress because a huge force is concentrated there in a very small area. Which may develop some permanent deformation at the contact surfaces. Therefore, this developed maximum stress must be less than a safe (or acceptable) limit equal to 4000 MPa. 14 Constraint 16: Rollers placed in between the rings may collide with each other while rolling. Therefore, to elude collision, an angular space of at least one degree is provided between the rollers.
Constraints 17&18: The roller effective length should be less than a particular value, B. Furthermore, it should not project out the bearing width and the chamfering of the rings. 17 Constraint 19: While rolling of the rollers and raceways, a shear stress called the maximum dynamic shear stress is developed below the surface of the raceways. Therefore, for a smooth running and to avoid failure of the bearing, the depth where maximum dynamic shear stress occur should be greater than the three times of the depth, where maximum static shear stress occurs.
The constraints from 1 to 19 that are discussed above must be converted into robust constraints using Eqs. (41) to (47) because they are still in the deterministic constraint format.
Constraints 20 to 25: This constraint is not applicable for deterministic optimization and applies only for the robust optimization. In physical programming for RDO, the designer is allowed to state the preference of deviation in objective functions, which may cause due to deviation in input design variables. As each design variable tolerance value input by us is the maximum possible tolerance, hence the objective function variation corresponding to all those design variable tolerance inputs collectively will also be a maximum possible, as can be seen in Eqs. (10) and (28). Also, the resulting variation in objective functions should always lie within the range prescribed by us and towards the lower limit since our aim is also to minimize the objective function variation. Therefore, constraints are While prescribing the ranges for percentage variations in the objective functions hmin with the help of this constraint, the percentage variation is prescribed for both hinner and houter. This will help us in ensuring that, while with the help of Eq. (28), we can reduce the magnitude of variation of any one among, hinner or houter, (depending on whether hmin is equal to hinner or houter ), we can also make the magnitudes of variation of both, hinner as well as houter, lie within a desirable range with the help of these constraints in the form of Eq. (47).

Variable bounds
To find the bounds in the case of deterministic optimization for basic design variables for the CRB NU 202, inputs which have been taken are as: boundary dimensions (from Table 1), the radial load equal to 10% of Cd (as mentioned in Table 1), the speed of shaft equal to the limiting speed for NU 202, i.e. 26,000 rpm (from Table 1) and the maximum contact stress (as mentioned in Table 2). As we know, the rated life (or L10) for a CRB is mathematically given as 10 3 10 where, L10 is in millions of revolutions, Cd is in N and radial load (Fr) is in N. Hence, we prefer to take a lower value of the radial load in order to achieve greater rated life for the CRB considering the fact that our aim in this paper is to maximize the fatigue life. Hence, we chose the radial load as 10% of Cd.
It must be noted that we have taken the speed of shaft equal to the limiting speed for each of the bearing, as it is unsafe to operate any bearing beyond its limiting speed. The bounds for the basic design variables can be derived from the constraints used in deterministic optimization. Hence, the bounds for Dm, Z, Dr and le for any CRB are  [25]. The bounds for these variables shall remain the same even in case of robust optimization as the case when they were being used for the deterministic optimization [25], because the variable bounds differ for any CRB for the robust and deterministic optimization only for those variable which have certain tolerances during manufacturing. Also, the bounds for these design variables would be the same irrespective of the type of CRB.
The bounds for the design variables of the CRB NU 202 have been mentioned for the deterministic optimization in Table 3. To convert the range of bounds for basic design parameters from the deterministic optimization to the robust optimization, we can use Eq. (40), but we will require tolerances for that. Tolerances are already available in the bearing catalogue and [2] has been referred for finding them out. Roller bearing engineers' committee (RBEC) Class 1 or normal class tolerances has been chosen as it provides maximum tolerances in the design variables. Defining acceptable tolerances for any dimension means that the bearing cannot have variations in dimensions beyond that value. The values of tolerances for the boundary as well as internal dimensions of each of the CRB have been mentioned in Table 4. While in the case of bearing boundary dimensions, we always take the tolerances as negative, here in case of basic design variables, we may take the tolerances as negative as well as positive. We take the input for tolerances in bearing boundary dimensions as negative because they have been conventionally accepted as negative values only. But it must be noted that the robust design constraints as well as the robust design objective functions treat negative as well as positive tolerances (whether of basic design variables or bearing boundary dimensions) in the same manner. Although the positive as well as negative tolerances are treated as same by the robust design constraints and objective functions, but we get different values of optimized design variables and optimized objective functions on taking a design variable tolerance as negative or positive. That is because of the fact that the sign of design variable tolerance will cause the new robust design bounds of basic design variables to be different even if, the tolerance magnitude remains same. Hence, although the final optimized objective functions and design variables for the cases of positive and negative tolerances shall have very close values, but the probability of them being exactly same is very less. However, the final optimized objective functions as well as design variables are not at all affected by the sign of tolerances in bearing boundary dimensions and only affected by their magnitude. Hence, the bounds for the design variables of the CRB NU 202 have been mentioned for the robust optimization for negative tolerances in Table 5.

Introduction to multi-objective deterministic optimization
Optimization is the process of achieving most suitable results for the maximization/minimization of an objective function with some associated restrictions. An engineering design optimization problem may be composed of one or more than one objectives, associated with a constrained parameter space. The optimization problem can be linear or non-linear, which depends on whether the objective functions and constraints are linear or non-linear, respectively. The present optimization problem of CRBs is non-linear.
While, the single-objective optimization problems (SOOP) consist of just one objective function, the multi-objective optimization problem (MOOP) refers to the solution of problems with two or more objective functions, which are at conflict with each other most of the time, i.e. increase in one objective leads to decrease in at least one of another objective or vice-versa. Thus, in MOOP, it is not possible to get a single solution that optimizes all the objectives at the same time. Fig. 2 represents the design parameter space and its corresponding objective space. For each solution x in the decision space, there exists a point in the objective space z. The mapping takes place between l-dimensional solution vector and m-dimensional objective vector. In MOOP, a vector of objectives is optimized, hence the MOOP is also known as the vector optimization. In MOOP, to compare two solutions concept of domination is used to check whether one solution dominates the other or not. A solution x is said to dominate (more optimal) another solution y, when the solution x is no worse than y in all objectives, and x is strictly better than y in at least one objective. Considering all the objectives are of the maximization type, the dominance of x is given as, where xy f means domination of x over y. The set of non-dominated solutions in a local search space is termed as non-dominated front. The global non-dominated solution set is the Pareto-optimal set. Such trade-off fronts are also termed as Pareto-optimal fronts (POF), named after Vilfredo Pareto who stated this concept in 1869. So, we can say that, all POFs are non-dominated but all non-dominated fronts are not POFs. The main task in MOOP is to get Pareto-optimal solutions (POS).
Evolutionary Algorithms (EA) mimics natural evolutionary principles to establish the search and optimization methodology. The typical feature of EAs is finding and maintaining multiple solutions in only one single simulation run. EAs deal with the population of solutions rather one solution in all iteration. The evolutionary algorithm has been chosen for this problem due to the above mention statement. Evolutionary algorithms available for MOOP are Non-Elitist MOEAs and Elitist MOEAs (Distance based Pareto GA (DPGA), Pareto Archived Evolutionary Strategy (PAES), Strength Pareto Evolutionary Algorithm (SPEA), SPEA-II, Non-dominated Sorting Genetic Algorithm (NSGA-II)). Out of these, NSGA-II has been chosen for this problem due to its low computational complexity and diversity preservation. The NSGA-II algorithm also exhibits elitism, but that can be seen in the most of other multi-objective optimization procedures as well.

Implementation of NSGA-II to CRB optimization
NSGA-II algorithm has been chosen as the optimization tool. The flowchart of its working procedure is given in Fig. 3. Nine design variables that are to be optimized here are encoded in real coded chromosomes and some modifications is done to the code to take the number of rolling elements as an integer value. Table 2 states the material properties of bearing, which are going to be useful as inputs for deciding the range for the basic design variables. Alloy steel 52100 grade has been chosen as the bearing material due to its high fatigue strength and wearing resistance, cost-effectiveness and long working life. After inserting the constraints and objective functions to the code, design variable bounds are inserted.
Computation is carried on i5-5200U CPU @2.20GHz Intel processor with 8Gb RAM using Dev C++. For all the runs, the population size is kept 1000. This large value is chosen to get more number of solutions in the POF. After finalising the population size, GA parameters (mutation probability, crossover probability, mutation index, crossover index and random seed) are chosen. The lower value of GA parameters gives the maximum spread of solution points and hence diverse results. Hence, we choose lower values for GA parameters [28]. The crossover and mutation distribution indices are used here to provide a set of widely spread solution points. Random seed is also used here to create initial population for each run, which are supposed to go further processes, like the crossover and mutation where new offspring is generated. List of GA parameters is given in Table 6. After finalising the population size and GA parameters, the code is run for a different prescribed variation in objective functions with a high value of number of generation. Generally, at such a large generation all the solutions get converged to their optimal value, so it is chosen high. The objective function variation for which we get the maximum value of Cd and hmin, is chosen as best variation for that bearing. Therefore, the variation in Cd and hmin for which there is maximum value of Cd and hmin is found to be 1.16-1.5 % and 0.1-0.2 %, respectively. The variations for all bearings are given in Table 7 that are to be a part of the Constraints 20 to 25. These variations are taken same for negative as well as positive tolerance because nature of tolerance does not affect the variation in objective function and that can be verified in Eqs. (10) and (28). It is noteworthy that this variation in objective function can be of any value as per the designer's requirement. Herein, it is chosen such a value which will give a maximum value of Cd and hmin.

Evaluate (P(t))
Assign rank using dominance depth method and diversity using crowding distance operator to (P(t))

t < T P (t) = selection(P(t)) (Crowded tournament selection) Stop
No Yes

Merge populations, P (t) = P(t) U P (t)
Assign rank dominance depth method and diversity using crowding distance operator to P (t) P(t + 1) = survivor(P (t) t = t + 1 Fig. 3. Working procedure of NSGA-II.

Results and observations
Robust optimization of Cd and hmin. is performed for NU 202 for prescribed variations of (1.16-1.5) % in Cd and (0.1-0.2) % in hmin. for negative tolerances for different number of generation and a population size of 1000. The resulting Pareto-optimal fronts (for the objective functions in our present case, i.e. f1 = {Cd -∆Cd} and f2 = {hmin -∆hmin} as we have to simultaneously maximize both {Cd -∆Cd} and {hmin -∆hmin}) for different number of generations are drawn simultaneously on a common graph, as can be seen in Fig. 4, where they are seen to be converging as the number of generations approaches towards 121. Hence, the final solution will be obtained from the Pareto-optimal front (POF) obtained corresponding to a population size of 1000 and a number of generations equal to 121, which has been shown in Fig. 5. The final solution (optimized dimensions) obtained from the POF for the robust optimization of Cd and hmin. for NU 202 for negative tolerances and prescribed variations of (1.16-1.5) % in Cd and (0.1-0.2) % in hmin. have been provided in Table 8. If the number of generations is increased further for the population size of 1000, there won't be much variation from the currently-obtained results, although the computational time would be increased. Convergence might not occur at the current value for number of generations, if the value for population size is varied.
Similarly, NSGA-II is implemented for other standard CRBs as well, for the negative as well as positive tolerances. The population size remains the same for each run (i.e. 1000), the number of generations is varied up to the point of convergence in design variables and objective function. The final optimized dimensions obtained for the robust optimization of Cd and hmin. for each of the CRB mentioned, for the negative as well as positive tolerances and prescribed suitable ranges of variation in the objective functions have been provided in Table 8.
For the final optimized dimensions obtained for NU 202 for negative tolerances, as given in Table 8, the bearing is drawn using a solid modelling software, where its cross-sectional view is shown in Fig. 6 and hence the geometrical feasibility is ensured and hence, it can be inferred that there is no interference between the dimensions.
As the optimized dynamic capacity from this work is found to be greater than the standard dynamic capacity available in the catalogue, indicating that bearing's life is improved. It can be Hence, the rate of change of Cd, hinner and houter with respect to Dm, Dr and le is positive. It means that for a positive tolerance in any of the dimensions Dm, Dr and le during the manufacturing of bearings, the objective functions will increase, while a negative tolerance in the same will cause a decrease in them. Also, the objective function Cd changes at the maximum rate for Dr and at the minimum rate for le. Although, a positive tolerance for the design variables Dm, Dr and le does not seem to be a problem, as it causes an increase in the objective function Cd, which is good for the bearing. But in that case, if we have not obtained the final optimized dimensions using the robust optimization and simply done the deterministic optimization, the tolerances may cause the final manufactured dimensions to be violating one or more than one of Constraints 1-19. Now, let us see the total percentage variation in all three, i.e. Cd, hinner and houter, due to the combined effect of tolerances in Dm, Dr and le. Hence, using results from Table 8 Therefore, from Eqs. (10) and (28) (99) Hence, as it was mentioned that the objective function variation for the design variable tolerance values as input given by us for any CRB, shall always be within the range prescribed by us in Table 7 and towards the lower limit. The same has been observed here, in Eq. (99) as well. From Table 7, the prescribed range of variation for CRB NU 202 is (1.16-1.5) % for Cd and (0.1-0.2) % for hmin. and we do have obtained herein a net variation of 1.16 % in Cd and 0.1966 % and 0.2 % in hinner and houter, respectively.
But, we must also check that the optimum results that we have got for the robust optimization is satisfying each of the constraints or not. For that, we have calculated the values for the left hand side of each of the robust optimization constraints corresponding to the results obtained for bearing NU 202 for negative tolerances. The obtained values for left hand side corresponding to each constraint has been mentioned in Table 9. Since all values are greater than zero, it can be clearly said that the optimum solution obtained for bearing NU 202 for negative tolerances and prescribed variation in the objective function, does satisfy all the constraints. Hence, we can say that, whatever results we are getting in Table 8, which are obtained based on robust optimization of Cd and hmin simultaneously for each of the CRBs, do satisfy each of the 25 robust optimization constraints. Table 9 Constraints violation for bearing NU 202 for negative tolerance.
Let us consider Constraint 1 and the value of its LHS in Table 9. The value of LHS for Constraint 1 is 10.5, which indicates that if we only have to follow Constraint 1 to find out the value of Dm then we can still decrease its value (currently obtained through optimization) by 10 hence, we get the optimum value for Dm, only after, each of the constraint containing Dm as a design variable, has been satisfied. It must be noted that each of the nine design variables, that we are considering in the present study may be present in a constraint as the only design variable, or it may be present as one of the several design variables in a constraint.
Similarly, if we consider Constraint 7, the value of its LHS is 1.7, which indicates that if we only have to follow Constraint 7 to determine optimal Dr and Kdmin., then we can still change the values of both of these (currently obtained through optimization) simultaneously or individually, such that, the value of the expression does not decrease more than 1.7 mm. But, since, in order to determine the value of optimum Dr and Kdmin we have to follow several other constraints in addition to Constraint 7; hence, we get the optimum value for each of the design variables, Dr and Kdmin, only after, each of the constraint containing Dr and Kdmin, as the design variable(s), in combined form or separately has been satisfied.
For any of the constraint (whether consisting of just one design variable or a set of design variables), the nearer is the value of LHS to zero, the more critical that constraint would be. This criticality of a constraint can be seen in reference to all the constraints or in reference to only those constraints, sharing at least one common design variable. If we consider the criticality of a constraint in reference to all the constraints, then, for the set of design variables present in this most critical constraint, the convergence to a stable value shall occur more quickly than convergence to a stable value for the set of design variables present in less critical constraints. The explanation for which has been given in the immediately next paragraph. If we consider the criticality of a constraint in reference to the constraints sharing at least one common design variable, then, out of all the constraints consisting of a design variable, the final optimized value for that design variable shall depend on the constraint that is most critical for that design variable.
In Table 9, for Constraints 11 and 20 through 25, the value of LHS corresponding to the optimum results for the robust optimization of bearing NU 202 for negative tolerances lies very close to 0, which indicates that these are the most critical constraints (in respect to all the constraints). Being a critical constraint means that for each run, these constraints will limit the range of design variables (out of the range that we have input) to optimize the objective function. Hence, the algorithm will keep only those populations as favourable parent populations in further generations, which do satisfy these constraints (i.e. the most critical ones). Hence, for the design variables that are present in these critical constraints, parent populations are carried forward for many generations and hence, for those design variables, the convergence to a stable value can be said to have occurred very quickly i.e., in the initial generations itself. Hence, in the present case of robust optimization of Cd and hmin of bearing NU 202 for negative tolerances and also for the robust optimization of each of the CRBs for positive/negative tolerances, design variables Dm, Dr, ε, Z and le, which are present in Constraints 11 and 20 through 25 shall converge to a stable value more quickly than the other design variables.
The final outcomes that we obtained for each CRB in Table 8, while some of the design variables like, le, Kdmin, Kdmax, ε, e and β are converging to the boundary values of their prescribed range, the others like, Dm, Dr, and Z converge to a value lying somewhere between the bounds of prescribed range. First of all, let us do the study of convergence of Dm, Dr, Z and le. As can be seen in Section 4, the upper and lower bounds of these variables are derived from the bounds of the same variables in the deterministic optimization, which themselves are based on a few of the constraints. As has been mentioned previously, if we consider the criticality of a constraint in reference to the constraints having at least one common design variable, then out of all the constraints consisting of a design variable, the final optimized value for that design variable shall depend on the constraint that is most critical for that design variable. Hence, for any of the basic design variables, convergence to a boundary value of the prescribed range shall occur only if, for that design variable, the critical constraint is same as the constraint that was used to find its upper and lower bounds for the deterministic optimization.
Coming to the convergence study of the remaining design variables. Consider Constraints 7 and 8 that consist Kdmin. and Kdmax and generate search space for roller diameter. From Eqs. (93) through (95), it can be observed that Cd and hmin have positive variation with respect to roller diameter (Dr), which concludes that with an increase in Dr the objectives Cd and hmin increases and vice-versa. Therefore, algorithm will try to find the maximum roller diameter from the constraints 7 and 8, and in order to choose that it demands a higher value of Kdmin. and Kdmax while keeping in mind that other 23 constraints are not violated. Thus the design variables Kdmin. and Kdmax are converging to upper limit of their specified range. Now discuss the convergence of design variable e that ensures the running mobility. From Eqs. (93) through (95), it can be observed that Cd and hmin have positive variation with respect to the pitch diameter (Dm), which concludes that with an increase in Dm the objectives Cd and hmin also increases and vice-versa. As Constraints 9 and 10 bound the search space for Dm through the design variable, e, and gives a maximum value of Dm when e has the maximum value. Thus e is converging to its upper limit, while keeping in mind that the remaining constraints are not violated. Now for the discussion of convergence of design variable, ε, consider Constraint 11, where a smaller value of ε increases the upper limit of search space of both Dm and Dr, to an even higher value. Hence, the optimization algorithm converges ε to the least possible value (i.e. lower bound for its prescribed range), meanwhile ensuring that the remaining constraints are not violated.
Again from Eqs. (93) through (95), it can be observed that Cd and hmin have positive variation with respect to the roller effective length (le), which concludes that with an increase in le, the objectives Cd and hmin increases and vice-versa. In Constraint 17, a higher value of β gives a higher bound of search space for le. Thus, β is converging to its upper bound, while keeping in mind that the remaining constraints are not violated.

Conclusions
In the present work, the multi-objective robust optimization is performed on a set of standard catalogue cylindrical roller bearings to give us an optimum design, which not only offers the maximum dynamic capacity and EHL film thickness, but also remains least affected by the possible practical error that can occur in the form of tolerances during its manufacturing. The study presented some mathematical formulations that make it easier to achieve our objective of having a robust design. The formulations have shown that in order to achieve the robustness, we need to make adjustments in all three components of any optimization process, i.e. the objective function, design variable bounds and constraints. Two objective namely dynamic capacity and EHL film thickness are optimized simultaneously subject to non-linear constraints using NSGA-II. The objectives are found to be conflicting and form a Pareto optimal front and the designer can decide which solution to take but the solution shown in the tables is the knee solution. Nonetheless, every solution on the POF is an optimal solution and one can choose any of the solution based on their preferences. Sensitivity analysis is also carried out to evaluate the objective function variation at the maximum possible tolerance in each design variable and results show that the variations are within their prescribed range. Important observations which are made from the results obtained are  A bearing is possible with such a design that its variation in performance due to tolerances in design variables (whose maximum value has been prescribed by us in positive as well as negative direction) during manufacturing, is always within the range prescribed by us. This is possible not only on paper but also can be manufactured, as solid modelling software can be used to draw such a bearing with geometrical parameters as obtained by the robust optimization.  The rate of change for both hinner , houter and Cd, corresponding to Dm, Dr and le undergoing variations during manufacturing, is positive.  Sensitivity analysis has shown that the objective function variation at the maximum possible tolerance in each design variable are within their prescribed range.
 The obtained POFs and POS are manifesting that objective functions are conflicting each other. All the solutions in the POFs are optimal. Therefore, designers can select any solution based on their requirement.  Out of all the parameters, which have been given tolerance, the roller diameter has the maximum effect while the roller length has the least effect on both the objective functions (Cd and hmin).
Hence, the present paper designs a rolling bearing that are robust in nature that means in addition to have maximum mean performance, it has least variation in its performance caused due to manufacturing errors. This design methodology can also be applicable for the advanced design of rolling bearings considering several design factors, like the radial internal clearance, crowning of rollers or raceways, heat generation, etc. Furthermore, it is not limited to the bearings only, but can also be used for the design of static or dynamic elements (like gears, levers, cams, etc.) of a machine. Load on the inner ring at the most heavily loaded roller (N) Qo

Nomenclature
Load on the outer ring at the most heavily loaded roller (N) r1 Outer ring chamfer height (mm) r2 Outer ring chamfer width (mm) r3 Inner ring chamfer height (mm) r4 Inner ring chamfer width (mm) , io e R Equivalent radius inner and outer respectively (m)

Availability of data and materials
All data and materials are in the manuscript itself.

Competing interests
There is no competing interests.

Funding
This research has not received any research funding.

Authors' contributions
The first author contribution is 50% (execution of whole formulation and coding and making draft version of the paper) and the second author contribution is 25% (making draft version of the paper, Writing -review & editing.) and third author contribution is 25% (In generation of idea, supervision of work, and fine tuning the paper)      Evaluate (P(t))

List of Figures
Assign rank using dominance depth method and diversity using crowding distance operator to (P(t)) t < T P (t) = selection(P(t)) (Crowded tournament selection) Stop

No
Yes Evaluate P (t) Merge populations, Assign rank dominance depth method and diversity using crowding distance operator to P (t)  Tables   Table 1 Parameters used for the CRB design consideration. Table 2 Bearing's material and lubricant properties. Table 3 Variable bounds for the deterministic optimization of bearing NU 202. Table 4 Tolerance (mm) in bearing's boundary and internal dimensions. Table 5 Variable bounds for the robust optimization of bearing NU 202 for negative tolerance. Table 6 List of GA parameters. Table 7 Prescribed ranges of variation in objective functions for each of the CRBs.   Evaluate (P(t))

List of
Assign rank using dominance depth method and diversity using crowding distance operator to (P(t)) t < T P (t) = selection(P(t)) (Crowded tournament selection) Stop No Yes P (t) = variation(P (t)) Evaluate P (t) Merge populations, P (t) = P(t) U P (t) Assign rank dominance depth method and diversity using crowding distance operator to P (t) P(t + 1) = survivor(P (t) t = t + 1