Rotor Design Optimization of Squirrel Cage Induction Motor - Part I: Problem Statement

Squirrel cage induction motor is the most widely adopted electrical machine in applications directly fed by the main grid. The analysis, design, and optimization of this machine topology has been addressed by a considerable amount of literature over the last century. Although its wide adoption, the induction motor design, especially when carried out in an automatic fashion, still presents significant challenges because the accurate prediction of the performance requires time-consuming finite element analysis. This work proposes a systematic approach to perform the design optimization of a squirrel cage induction motor focusing on the rotor slot geometry, being this the major player in defining the torque-speed characteristic. Structured as a two-parts companion papers, this first part presents an innovative performance evaluation methodology which allows a very fast estimation of the torque, and efficiency behaviour preserving the results’ accuracy. The proposed performance estimation technique is assessed against experimental tests carried out on an off-the-shelf induction motor. The selection of the performance indexes to be optimized is justified in detail along with the description of a generalized rotor parametrization which allows a comprehensive exploration of the research space. The presented optimization procedure is then applied to a case study, and the preliminary results are commented, highlighting benefits, and drawbacks of the proposed methodology.


I. INTRODUCTION
S QUIRREL cage induction motor (SCIM) is the most widespread electrical machine for historical reasons, manufacturing simplicity, robustness and low cost. Typical industrial applications for directly fed SCIMs include fans, pumps, compressors, conveyors, presses, elevators, extractors, etc, with the first three covering the majority of the market [1]. National and international standards provide a set of minimum performance requirements for general purpose motors. The National Electrical Manufacturers Association in US (NEMA) categorises SCIMs, having the same power, in different typologies [2] Manuscript received June 29, 2020; accepted July 29, 2020. Date of publication August 27, 2020; date of current version May 21, 2021   according to the maximum value of the starting current and the minimum values of the starting, pull-up and breakdown torque (as shown in Fig. 1(a)). An overview of the NEMA classification of small-medium SCIM (power rating: 0.367-367 kW) is reported in Table I, where torques and currents are expressed in percentage of the respective rated values. In Fig. 1(b), a qualitative representation of the torque-slip curves of each class is reported. Design classes A and B exhibit a similar trend in terms of performance requirements and they are used whenever the starting torque demanded by the application is relatively low. Classes A and B differ from each other in terms of the starting current requirement: the first one does not define any limit and it is employed for higher efficiency and breakdown torque applications with respect to the design B. Such motors (class A) typically require the use of a reduced starting voltage strategy. They are typically designed with simple rotor slot geometries while class B motors are realized with double or deep rotor slot configurations. Design class C imposes a higher starting torque requirement than design class B with the same limitation in term of starting current, at the cost of a lower breakdown and pull-up torque and a higher rated slip and a lower rated efficiency. They are typically realized with double rotor slot configuration with rotor resistance higher than motors design class B. It is used whenever starting capability under high load is required. Motor's design of class D, usually realized with single rotor slot geometries, guarantee a really high starting torque at the cost of a very high rated slip and low rated efficiency. They are used in applications with high peak loads, high inertia and very intermittent operation [3]. Table I reports the range of variation of each requirement divided in classes; while in Fig. 2 is reported the minimum starting, breakdown and pull-up torques as a function of the power rating of the motor, along with the maximum starting current for a 4-poles SCIM class B (which covers about 90% of all general purpose induction motors [4]).
Similar classification and minimum performance requirements in terms of torque-speed characteristics are defined by the IEC 60034-12 however, for the sake of brevity, are not reported in this paper [5].
Before the energy crises in the 1970s, the majority of SCIMs for general-purpose applications were designed to provide the rated power and operating characteristics at the minimum cost; the efficiency was not a strict requirement, but rather a consequence of the maximum temperature limit. As the energy price began rising, motor manufacturers started promoting motors with higher efficiency. Regulatory authorities begun introducing minimum efficiency requirements, providing classification of electric motors also according to this figure of merit [6].
The NEMA defines two categories of efficiency: the energy efficient and the premium efficiency. In Fig. 2, the minimum requirement for both the classifications is reported. Similarly, the IEC 60034-30-1 defines four different classes of efficiency (they are not reported for the sake of brevity), where two of them are essentially equivalent to NEMA levels: IE2 (Energy Efficient) and IE3 (Premium Efficiency).
In Fig. 2, the gaps between efficiency classes are quite narrows, therefore the evaluation of the SCIM performance during the design process requires high accuracy to appreciate the improvements. The most accurate evaluation of the SCIM performance can be achieved with computational expensive time step finite element analysis (FEA). However, a very precise and time consuming FEAs is not suitable when an iterative design procedure is adopted [7]. On the other hand, a fast performance prediction based on analytical formulation might lack in accuracy due to the inevitable assumptions of the latter. Although both stator and rotor geometries affect the SCIM performance, the shape of the slot hosting the rotor bars, plays the most important role in determining the behaviour of the torque-speed characteristic of the machine as well as its efficiency as a function of the load. Traditionally, the rotor slot design process starts from a limited set of well-known slot shapes (trapezoidal, rectangular, rounded, Boucherot) and geometrical modification are then applied until the considered performance specifications are met [8], [9]. By doing so, the quality of the obtained design largely depends on the engineer experience and on the tools used to predict the SCIM performance.
An attempt to ensure that the obtained design is not a local optimum, given the performance indexes and constraints, is to use a formal optimization algorithm. Several SCIMs' optimization studies have been proposed in literature. Some of them [8], [10] employ deterministic optimization algorithm using aggregated cost function while others [3], [11], [12] utilize stochastic search methods targeting specific applications. All the proposed works lacks in generality, either due to the selected rotor slot parametrization or the chosen objective function or the adopted performance evaluation method (which sometimes not include all the involved non-linearities or it is computationally too expensive to allow a complete exploration of the research space).
This work, structured as two-parts companion papers, proposes an automatic design methodology for SCIMs. Although the procedure is general and potentially suitable to design the whole machine, in these papers it is applied to optimize only the rotor slot geometry, being the most critical part in defining the torque and efficiency performance.
In this first part of the work the optimization problem is defined as follows. First, a detailed description of an innovative mixed analytical-FE performance evaluation method is presented in Section II, along with the experimental validation in Section III. It will be demonstrated that the selected approach guarantees obtaining a good trade-off between the conflicting requirements of high accuracy and short computational time. Afterwards, Section IV describes the optimization process and the selection of the objectives and constrains. Finally the iterative procedure implemented to relieve the computational burden of the performance identification is presented. A generalized rotor slot parametrization, enabling to explore a wide variety of rotor slot shapes during the optimization, is then described in details. Section V and VI report the preliminary optimization results and the final considerations highlighting the benefits of the proposed optimization procedure.
The second part of the work will outline a detailed description of the various optimization results. Different optimum designs obtained considering several constraints and rotor cage materials will be commented considering also the thermal behaviour of the machine.

II. PERFORMANCE EVALUATION METHOD
In general, the performance calculation methods of a SCIM can be divided in two main families. The first one relies on brute force time-step finite element analysis (FEA), while the second is based on the resolution of the single-phase equivalent circuit (EC) [7], [13], [14], reported in Fig. 4. The first method tends to be more accurate but computationally expensive, while the results quality of the EC-based methods depends on the accuracy in estimating the unknown lumped parameters. The latter can be calculated analytically [15], [16] disregarding the parameters' dependency with respect to current and frequency or alternatively from a set of FEAs. Several procedures have been proposed in literature to estimate the EC parameters as a function of current and frequency with less expensive magneto-static or time-harmonic (TH) FE simulations [7].
The most common method consists in FE-replicating the indirect experimental tests, i.e. no-load and locked-rotor tests. In particular, the no-load test, carried out by magneto-static or time-harmonic FEAs for different current values, provides the equivalent no-load inductance and iron losses resistance. The subdivision between magnetizing and leakage components of the inductance can be performed through the calculation of the first harmonic of the airgap flux density or via a second set of simulations, supplying properly the rotor bars [17]. The locked-rotor tests, executed via standstill TH simulations for different value of the frequency in the rotor reference frame, aims at determining the rotor impedance. The equivalent impedance resulting from this test can be determined either from the FE-calculated rotor joule losses and the total magnetic energy or applying Kirchhoff's circuit laws with the FE-calculated voltage and current phasors [18]. The rotor impedance is then obtained from simple circuital considerations knowing the magnetizing and the stator leakage inductance. This classical approach presents two critical issues, strictly correlated with its implicit assumptions: r the EC parameters determined with the FE replicated indirect tests do not change in the real load scenario; r the parameters calculated from the no-load tests depends only from the current while the one determined from the locked-rotor test depends only from the frequency; this is not true in the case of closed rotor slot geometries. Although the last point can be overcome by identifying the EC parameters with FEAs as function of both current and frequency, at the cost of an increased computational burden (performing non-linear TH-FEAs), the first issue is an intrinsic limitation of this procedure.
Recently, a different method to calculate the EC parameters has been presented [17], [19] with the aim of solving both limitations. This alternative approach employs one standstill (non-linear) TH-FE simulation, carried out for a single load condition, i.e. imposing a certain current and frequency in the rotor reference frame, to determine all the EC parameters related with such operating point. By doing so, all the EC parameters corresponding to one operating point are simultaneously determined from one FEA, replicating as close as possible the real working scenario. This method, named contextual identification, relies on the knowledge of both amplitude and phase of the equivalent transformer ratioK, which can be determined analytically or from a single no-load TH-FE simulation [17].
For each operating point, simulated with a current fed standstill TH-FEA, the magnetizing currentİ m is the sum of the imposed stator currentİ s and the FE calculated current flowing in the reference rotor barİ r divided by the complex ratioK. The stator leakage L ls and magnetizing M inductance can be then determined as: whereψ s is the FE-calculated stator flux phasor, while the equivalent iron loss resistance is given by: where P fe are the iron losses resulting from the TH-FEA (the stator part of the iron losses is multiplied by 1/s). Once the parameters on the primary part of circuit have been calculated, knowing the complex ratioK, the voltage induced on the reference rotor bar can be calculated as follows: Finally, from the secondary voltage phasorV mr , the rotor parameters can be easily computed as in (5) and (6): Fig. 3 reports the EC parameters calculated with the above identification method as function of both current and frequency for a SCIM featuring closed rotor slot; clearly all four parameters depends on both current and frequency. Once such parameters are identified in the current-frequency plane, and the 3D parameters (L ew , L ring , R ring ) have been analytically calculated [15], [16], the non-linear EC can be solved for any slip value in order to fully determine the steady state performance of the SCIM. Although the contextual identification method allows a more accurate estimation of the EC parameters of any SCIM, including designs featuring closed rotor slot geometries, it is computationally expensive and therefore not suitable to be embedded in an automatic design procedure. The computational time is directly proportional to the number of currents chosen to identify the behaviour of the EC parameters for a given frequency, i.e. slip. The choice of the interpolating functions, used to approximate the parameters in the resolution of the EC, is important because it can lead to a reduction of the number of currents to be simulated. However the selection of these interpolating function is crucial and becomes critical when very similar machine geometries are compared, as in the case of automatic design procedure. With the twofold aim of reducing the computational burden while ensuring high accuracy during the EC parameters calculation, an iterative procedure is proposed in this work. The latter, shown in Fig. 5, for a given slip, consists of: r carrying out the time-harmonic FE simulation (with a supposed stator current I s,0 at the first iteration); r identifying the EC parameters via the Eq. (1)-(6); r solving the non-linear EC and so determine the effective stator current I s,k+1 . If the relative error between the FE simulated stator current I s,k and the one calculated solving the EC I s,k+1 is below   a certain limit I , the iterative procedure ends, otherwise a new FEA is carried out with the new value of the stator current I s,k+1 . This procedure permits to reduce the number of currents to be simulated at a certain operating point and to avoid the use of interpolating functions required during the EC resolution.
Implementing this iterative method permits to reduce the number of simulated currents per slip and the overall computational burden, ensuring high accuracy of the evaluated performance.

III. EXPERIMENTAL VALIDATION OF THE PERFORMANCE EVALUATION METHOD
For validation purpose, the performance of an off-the-shelf SCIM are computed by means of the proposed methodology and compared against experimental results. The main characteristics and geometrical parameters of the motor are reported in Table II. The validation exercise aim to asses the accuracy of the approach in predicting the steady state behaviour of the machine. Fig. 6 shows the experimental test rig set up, consisting in: the SCIM under test on the right side of the picture connected through the torque transducer to a DC motor acting as a prime mover. The DC motor operate in speed control while the SCIM is fed with nominal voltage and frequency. Once the steady state temperature of the motor at rated operating point is achieved, different slip values are investigated in order to identify the steady state performance of the motor under test.
Since the SCIM features a skewed rotor layout, a modification of the EC is required. The approach presented in [20] is adopted  and the modified circuit is reported in Fig. 7 for the sake of clarity. The coefficient K sk compute as in (7) is used to weight the magnetizing inductance and to split part of it on the stator and rotor branches of the circuit: where θ sk is the electrical skewing angle. The EC lumped parameters are computed using 2D FE model with straight rotor cage. Some of the parameters are then modified according to Fig. 7 and finally the performance in term of steady state torque and stator current are computed by solving the non-linear EC as described in Section II. Fig. 8 shows a comparison between the experimental and the predicted results in per unit. The compared quantities are showing a good agreement, thus confirming the accuracy of the adopted method.

IV. OPTIMIZATION PROBLEM DEFINITION
The following three sub-paragraphs detail A) the selection of the objectives and constrains of the optimization problem, B) the iterative procedure implemented to relieve the computational cost of their identification, and C) the rotor slot parametrization based on a generalized Boucherot bar.

A. Objectives and Constraints
Starting torque and current, pull-up and breakdown torques, rated efficiency and power factor can be all considered possible objectives or constraints of the optimization problem under study. In fact, all of them are deeply influenced by the rotor slot geometry and their values have to comply with the boundaries defined by the national or international standard (except the power factor). If all the above mentioned performance indexes are determined, the subdivision between objectives and constraints depends on the specific application under study. For example, applications where the rated efficiency is the most important requirement, the latter can be treated as the only objective of the optimization problem while all the remaining indexes can be considered as constraints. However, the accurate evaluation of all of them is computationally expensive because it implies the determination of the full torque-speed characteristic.
In this work, with the aim of reducing the computational burden and the complexity of the problem, only the rated efficiency, the starting torque and current are determined during the optimization process. The first two (η rated , T start ) are selected as objectives while the ratio between the starting and the rated current (k I ) is considered as a constraint; the optimization problem is therefore defined as: Selecting the current ratio as a constraint and not as the third objective reduces the complexity of the optimization problem and so increase the quality of the results [21].

B. Evaluation of the Starting and Rated Performance
The starting performance are determined carrying out a single voltage-fed TH-FEA. In fact, in such operating point, the iterative procedure based on current-fed time harmonic FE simulation, described in Fig. 5, is not needed. At starting condition, s = 1, the rotor is at standstill and so the pulsations of stator and rotor quantities are equal. Consequently, a single TH-FE simulation supplying the stator windings with symmetric threephase voltages, with rated amplitude and frequency, allows the direct evaluation of the performance. The rated operating point is hereafter defined as the condition in which the motor provides the rated torque T rated . However, the rated slip s rated is not a-priori known as it strongly depends on the rotor slot geometry. To identify the effective operating point, an iterative procedure, summarized in Fig. 9, has been implemented, consisting in the following steps.
i) The values of the rated slip, rated current and the slopes of both current and torque function at the rated condition (s, I s0 , m T , m Is ) are initially supposed or determined following the procedure described at the end of this section. ii) The machine performance are evaluated using the iterative procedure reported in Fig. 5 hereafter called TH-FE current iteration. iii) The procedure ends if the obtained torque T (i) is equal to the rated torque within the limit defined by T such as |T (i) − T rated | < T . iv) Otherwise, the slip and current values are updated according to (9) and (10) using the results of the last simulation and the slopes m T and m I s : v) The machine performance are reevaluated by means of the TH-FE current iteration. vi) If |T (i) − T target | < T , the rated operating point is identified. vii) Otherwise the slope of both current and torque are updated by means of (11) and (12): and the new current and slip are determined according to (9) and (10). viii) The machine performance are reevaluated again.
ix) The iteration stops if the inequality |T (i) − T target | < T is satisfied, otherwise the procedure restart from step vii.
At the end of every iteration the rated performance (s, I s , m T , m Is ) corresponding to a certain rotor geometry, described by the vector x, are stored. During the optimization, when evaluating the performance of the geometry x * , with the aim of reducing the number of FE simulations carried out during the search of the rated current and slip, the initial values of the variables used to start such iterative procedure are imposed equal to the corresponding variables of the closest geometry previously evaluated x p min . In other words, before starting the evaluation of the geometry x * , the euclidean distances between the latter and all the previously evaluated solutions x p are determined: where n v is the number of variables defining the rotor geometry. Once all these distances are calculated, the nearest solution x p min to the current one, can be easily identified as well as its rated performance (which have been opportunely stored). In the next section, it will be shown that adopting this "learning" ploy the number of FE simulations needed to identify the rated performance of a certain geometry is greatly reduced. This is because the evaluation starts from the current and slip values of the most similar geometry previously evaluated.

C. Rotor Bar Parametrization
A generalized parametrization of the Boucherot bar has been implemented in order to widen the shape variety explored during the optimization. The rotor slot geometry is described by 7 variables, as shown in Fig. 10, i.e. the neck width d, the radii of the outer, central and inner circles R out , R m and R in , the distance between the outer and central circles h out and the central and the inner circles h in and the external bridge thickness b. All these 7 variables are defined in per unit of their respective maximum values. In particular the following drawing sequence has been adopted: i) the bridge thickness b is calculated knowing its per unit value b pu and its upper b max and lower limits b min : ii) the outer circle radius can be defined once its maximum value has been calculated from geometrical consideration: (15) iii) the distance between the outer and central circles follows after the calculation of its maximum allowed value: iv) while the distance between the central and inner circle is: v) once the center of the central x m and inner circles x in are defined, the maximum value of their radius are calculated  imposing geometrical constraints: thus allowing to calculate the radii R m , R in from their per unit values; vi) finally the neck width can be calculated as: Adopting this per unit parametrization allows obtaining always geometrically feasible solutions, and to explore a wide variety of bar geometries including Boucherot types, drop shapes and deep bars as shown in Fig. 11 . Varying without any constraint all the 7 variables describing the geometry implies that the slot surface A slot is not constant. If the latter has to be kept constant (A slot = A target ), the variables are not anymore independent and a correlation among them has to be established. Given the difficulty of finding a closed form of such relationship, the variables defining the slot geometry (all except b) are updated proportionally to the square root of the ratio between the target slot area and the actual one A target /A slot . By doing so, the area of the slot equals the target value and its shape is preserved respect to the initial geometry.

V. PRELIMINARY OPTIMIZATIONS RESULTS
A multi-objective stochastic optimization algorithm (NSGA-II embedded in Matlab) has been chosen to carry out the rotor slot design of the SCIM with a population size of 100 elements  evolving for a maximum of 100 generations. As previously mentioned, the stator geometry has been kept constant during the optimization process and its main geometrical parameters are reported in Table. II. Figure 12 reports the results in terms of Pareto fronts of two optimizations carried out considering a maximum current ratio k I,max of 7.4 and a constant slot area A * slot . The first optimization considers the full set of geometrical variables during the design process (called variable bridge in Fig. 12), while in the second one the bridge thickness is kept constant and equal to the minimum allowed by manufacturing considerations (0.5 mm). The analysis of the results is leading to the following considerations: r rated efficiency and starting torque show a clear competitive behaviour, i.e. the improvement of one of them implies the worsening of the other; r the range of variability of the starting torque is much higher than the respective range of the rated efficiency; e.g. an increment of 100% of the starting torque (from 140 to 280 Nm) implies a reduction of the rated efficiency of only 1.1%; r both optimizations lead to almost the same results in terms of rated efficiency and starting performance; r the real current ratio of the optimal solutions is substantially equal to the maximum value allowed during the optimization process. Analysing means (μ) and standard deviations (σ) of the geometrical variables of the optimal solutions, obtained considering the whole set of parameters during the optimization, shown in Table III, it is evident that the bridge thickness tends to converge to the minimum value while the distributions of other variables have a much higher variance. The number of geometrical variables to optimize can then be reduced, keeping the bridge thickness equal to the minimum value, simplifying the optimization problem without affecting the results (as shown in Fig. 12).
The Pareto fronts have been obtained evaluating a total number of 10 000 geometries. Each solution requires the assessment of the starting performance via a single voltage-fed TH-FEA and a variable number of current-fed TH-FEAs to determine the rated performance. Fig. 13 reports the evolution of the average (per generation) number of the simulated slips and the total TH simulations performed during the search of the rated operating point. It can be observed that the number of simulations steeply decrease at the fifth generation when the learning ploy, described in Section IV, is activated. Once enabled, the iterative procedure to identify the rated performance of the solution under assessment starts from the current and slip values of the most similar geometry previously evaluated. The number of simulations continues to reduce as the optimization algorithm proceeds (as shown in the inset of the same figure) because more and more results are available to be compared with the geometry under evaluation and so the probability of having similar slots increases. After few generations, each solution requires on average 1.25 slips to be simulated in order to find the rated condition. Each of them requires on average less than 1.5 TH-FE simulations to identify the rated current.

VI. CONCLUSION
This work, first part of two companion papers, has proposed an automatic design procedure of squirrel cage induction motors including the following novelties.
r First an original performance evaluation method has been presented. The latter is based on an iterative resolution of the classical single-phase equivalent circuit featuring parameters determined via current-fed finite element time harmonic simulations in the real load scenario. The presented method allows a really fast prediction of both torque and efficiency without loosing accuracy as verified by the experimental tests carried out on an off-the-shelf induction motor.
r Albeit the adopted performance evaluation method is general and so suitable to analyse and design the whole machine in a very time-efficient manner, as a vessel to investigate its capabilities, it has been applied to optimize the rotor slot geometry. A generalized parametrization of a Boucherot bar has been then introduced which allows exploring the most common slot shapes (double cages, drop-like, deep bar, etc.) with the minimum number of geometrical variables.
r Starting and rated performance in terms of torque and efficiency are the selected objectives to be optimized during the proposed automatic design procedure, having considered only applications directly fed by the main. Being the rated operating point dependent from the slot geometry, a new iterative algorithm has been implemented for its identification. The latter, essentially based on the secant method, starts the identification using the rated values of both slip and current of the most similar geometry previously evaluated. Adopting this approach, the number of FE simulations needed to identify the rated performance of the machine under assessment is greatly reduced. It has also been shown that this computational burden continues to decreases as the optimization algorithm proceeds because the likelihood of finding similar geometries increases. Thanks to this method, the evaluation of the objectives for each induction motor requires 1 TH-FEA for the starting performance calculation and about 1.3 TH-FEAs for the identification of the rated condition. Considering that each TH-FEA takes approximately 4 s on a standard-spec PC (Intel Core i3, 3.4 GHz, 8 GB RAM), the whole optimization process, which evaluate 10 000 geometries, takes less than 32 hours.
Part II of these companion papers will analyse: r several optimization results obtained considering different constraints in terms of maximum starting current, rotor slot area and cage materials; r the accuracy of the efficiency evaluation method implemented during the optimization, which does not consider the influence of the spatial field harmonics in the losses calculation; r the thermal behaviour of the optimized machines and its influence on the optimization results, which have been obtained considering a constant temperature distribution.