Improved Design-Oriented Analytical Modelling of Switched Reluctance Machines Based on Fröhlich-Kennelly Equations

The need to reduce computational burden in the preliminary design stage of switched reluctance machines is fostering the interest in design-oriented analytical modelling. To this end, this work proposes a novel design-oriented analytical model that comprises two main parts, each of them containing an original scientific contribution: 1) a new interpolation technique for the flux loci based on 2${\text{nd}}$-order Fröhlich-Kennelly equations, and 2) an analytical model that calculates flux linkage in partial overlap and saturated conditions. The model is validated against finite element analyses of four switched reluctance machines and the experimental results of a physical prototype. Finally, an in-depth discussion on the use, accuracy and limitations of proposed analytical model is provided.

Roberto Rocca was with the Department of Astronautical, Electrical and Energy Engineering (DIAEE), Sapienza University of Rome, 00164 Rome, Italy.He is now with the Department of Electrical Systems, CIRCE Technology Centre, 50018 Zaragoza, Spain, and also with the CIRCE Mixed Research Institute, (CIRCE Technology Centre and University of Zaragoza), 50018 Zaragoza, Spain (e-mail: rrocca@fcirce.es).
Paolo Bolognesi is with the Electrical Machines, Power Electronics and Drives Group, DESTEC Department, University of Pisa, 56122 Pisa, Italy (e-mail: paolo.bolognesi@ing.unipi.it).
Color versions of one or more figures in this article are available at https://doi.org/10.1109/TEC.2023.3307574.
Digital Object Identifier 10.1109/TEC.2023.3307574where a lower accuracy may be accepted in return for a faster execution [9].In a bid to provide a reliable alternative to FEA, designoriented analytical models should be developed to ensure a high accuracy for SRMs featuring any kind of geometry and conceived for different rated operating conditions.The analytical design of any SRM starts from the definition of its flux linkage loci, i.e., a set of curves expressing the phase flux linkage Ψ ph as a function of both phase current i ph and rotor position θ (Ψ ph (i ph , θ)), from which the rest of the performance is calculated [10].A schematic example is shown in Fig. 1.The process to define it analytically is usually divided into two sub-problems: 1) Choice of an analytical expression of Ψ ph (i ph , θ), e.g., Polynomials, Fourier Series, etc. 2) Analytical evaluation of the flux linkage values at the desired nodes (pair of values of i ph and θ).In general terms, the main challenge underpinning both points above is modelling the nonlinear effects of saturation.With regards of the analytical expression of Ψ ph (i ph , θ), two approaches are commonly followed, namely Fourier series and piece-wise expressions.For what is concerned with Fourier series, many contributions agree that this approach suites well the development of control systems, whereas its complexity is not ideal for design purposes [11], [12], [13], [14], [14], [15], [16], [17].In terms of piece-wise expressions, a relatively simple approach is to interpolate Ψ ph (i ph , θ) by two-linear-piece ones [18], with the issue of low accuracy nearby the knee-point of the saturation curves.Alternatively, a three-nonlinear-piece expression is proposed by T.J. Miller in [9], where a set of Ψ ph vs θ curves for constant i ph is defined (see Fig. 1), whose nonlinear parts are interpolated based on 1st-order Fröhlich-Kennelly (FK) equations [19].Although this approach has already proved to suite well the needs of designers, [20], [21], it still leaves room for improvement.Indeed, when interpolating a nonlinear region (e.g., between θ 2 and θ u in Fig. 1), a 1st-order FK equation only offers three degrees of freedom, resulting in non-smooth curves.Moreover, in [9], all rotor positions ranging from θ a to θ m (see Fig. 1) belong to the same "piece", posing the risk of losing accuracy in the vicinity of θ 1 .In bid to tackle the two above issues, this work proposes an innovative piece-wise formulation based on 4 pieces and whose nonlinear ones are expressed through a 2nd-order FK equation.
Moving on to the analytical evaluation of the flux linkage at the desired nodes, the challenge lies in modelling the saturation effects in partial overlap conditions.In [22], this situation is accounted for by dividing stator and rotor teeth into a saturable tip zone and a non-saturable trunk zone, although no relationships are defined between these zones, the rotor position and the saturation level.Alternatively, [23], [24], divide the stator tooth into a saturable and a non-saturable straight segment, whose thickness changes with the rotor position.However, despite its good accuracy, this approach results quite complex for design purposes.In [20], [25], [26], a multi-branch magnetic circuit is proposed, where flux tubes are shaped by tracing contour plots obtained via FEA, resulting in a considerable first-implementation effort.Moreover, flux tubes shapes may change significantly from one SRM to another, undermining the versatility needed in a design process.As an alternative, [27] develops a 3-reluctance network with one representing the teeth saturation.However, also in this work no relation between geometry, rotor position and saturation level is defined.The same inconvenience is found in [28], [29], where a basis and a winding function approach are respectively proposed.
In a bid to fill these gaps found in the literature, this work proposes a novel design-oriented analytical model catering for the Ψ ph (i ph , θ) loci, which is based on a twofold scientific contribution: 1) An improved interpolation technique of the Ψ ph (i ph , θ) loci based on a four-piece expression, whose nonlinear ones are expressed through a 2nd-order FK equation, 2) A novel analytical model catering for the flux linkage in saturated cores in partial overlap conditions.This article is organised as follows.Section II provides an overview of design-oriented analytical modelling of SRMs.Section III presents the unsaturated inductance profile.Sections IV and V are devoted respectively to the description of the FK-based formulation of the Ψ ph (i ph , θ) curves and the flux linkage evaluation in saturated cores and partial overlap.Finally, Section VI validates the proposed model, while Section VII provides a final discussion.

A. General Background
An example of SRM candidate is provided in Fig. 2, whose geometric variables are listed in Table I.The first step in modelling is usually the definition of the SRM Ψ ph (i ph , θ) loci.Subsequently, under the widely accepted assumption of mutually decoupled phases, the static torque loci of one phase T sta (i ph , θ) are evaluated by means of the energy/co-energy theory as in (1) where the integral refers to one stroke undergone by one phase [30]: Then, by keeping the assumption of magnetically decoupled phases, the phase voltage equation can be solved to obtain the variation over time of i ph and θ during one stroke [10]: Here, v ph denotes the instantaneous phase voltage imposed by the converter, R ph the phase resistance and ω the speed of rotation, which is normally kept constant at a desired design point.Then, as (2) is solved, i ph (t) and θ(t) can be plugged into (1) to obtain the instantaneous phase torque waveform T ph (t).Finally, contributions provided by the remaining phases can be found by properly shifting i ph (t) and T ph (t) with respect of θ.
Regarding the hypothesis of negligible mutual coupling effects, it is reminded that this last is normally accepted in the preliminary design stage, although its verification is necessary at later stages.

B. Design Candidate Performance Prediction
In the case of this work, the background described above is implemented in the design process shown in Fig. 3.As it can be seen, once the SRM candidate geometry is generated (this aspect is not addressed in this work), the first step to take is the phase unsaturated inductance vs. rotor position profile evaluation.This task is described in Section III.Once the unsaturated inductance profile is determined, the flux linkage vs. phase current at fully aligned and mid-way positions are calculated (see curves marked in red in Fig. 4) in the way outlined in Section IV.As it can be seen in Fig. 3, aligned and mid-way saturated flux linkages feed the "Flux Loci Calculator".Here, the flux linkage loci are evaluated through the interpolation process described in Section V.At this point, as Ψ ph (i ph , θ) loci are available, both (1) and (2) are resolved in the "Performance Calculator" block, so that the RMS current i ph−rms , average and torque ripple T AV G , T RP L are obtained.

C. First-Order FK Equation
In order to interpolate nonlinear SRM characteristics, e.g., between θ a -θ m or θ 2 -θ u in Fig. 1, [9] proposed the use of 1storder FK equations.Their generic mathematical formulation, A(θ), is as follows: (3) As it can be seen, A(θ) possesses four degrees of freedom, namely a x , b x , θ x and A 0 , which need to be set in accordance with the boundary conditions shown later in Section III.The key observation concerned with (3) is that a boundary condition setting a derivative equal to zero, as per maximum alignment or misalignment conditions, would force either a x or b x to be zero, which would be a non-sense (a x is in fact a combination of the remaining three parameters [21]).This, in turn, causes A(θ) to be non-smooth, resulting so in inaccurate profiles representations as shown later in Section V-E.To overcome this issue, this work proposes a 2nd-order FK equation, which is described in Section III-C.

A. Profile's Regions Definition
An example of profile of incremental self-inductance of one phase vs. rotor position L u ph (θ) is shown in Fig. 5.The idea proposed in this work is to model this curve by dividing it into four regions.The initial rotor position θ a is considered at the full alignment between one stator and one rotor tooth, as shown in Fig. 6(a).As the rotor moves, two facing teeth remain fully overlapped until the tip-to-edge position θ 1 (Fig. 6(b)), where Region I terminates.Then, Region II begins as stator and rotor partially overlap.The centre of the partial overlapping zone is referred to as mid-way position θ m , where just half stator tooth Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.overlaps with the rotor (Fig. 6(c)).From θ m onward, the rest of the partial overlap zone defines Region III, which ends at the tip-to-tip position θ 2 (Fig. 6(d)).Finally, Region IV describes the non-overlap zone up to the full unalignment condition θ u , which is sketched in Fig. 6(e).

B. Unsaturated Inductance Profile
To define the unsaturated inductance profile, the initial task is the evaluation of the unsaturated inductances corresponding to the rotor positions mentioned above: In this work, L u a , L u 2 and L u u are evaluated from the SRM geometry through the analytical technique proposed in [30], which is based on a permeance approach representing the magnetic flux tubes with an elliptic shape.Besides, L u 1 and L u m are calculated by the well-established assumption that in absence of saturation the inductance profile is constant between θ a and θ 1 and descends linearly from θ 1 to θ 2 : At this point, all regions can be defined.Thanks to the assumption of neglecting fringing and rounding phenomena, Region I is modelled through a constant inductance value equal to L u a .Subsequently, Regions II and III are modelled together assuming that the profile decays linearly until position θ 2 .Finally, Region IV begins, whose profile is obtained via the 2nd-order FK equation described in the next subsection.

C. Second-Order Fröhlich-Kennelly Equation
As mentioned before, for modelling the nonlinear regions, this work proposes a 2nd-order FK equation, whose expression is shown below: As opposed to its 1st-order counterpart, this formulation possesses four degrees of freedom: L x , a x , b x and θ x , which, in turn, allow one to set four boundary conditions, namely the values of the inductance at both region extremes, along with the derivative values.In the case of Region IV in unsaturated conditions, boundary conditions in (7) apply.
Thanks to this formulation, the first derivative constraint can be set that imposes a zero derivative at full unalignment, while the second one sets the slope at θ 2 equal to that of Region III.Once the boundary conditions are defined, a system of four equations can be written and solved for L x , a x , b x and θ x .
It is finally remarked that 2nd-order FK equations are used in this work for modelling all nonlinear regions, including in saturated conditions.This aspect is discussed later in Section V, where also the need to define the unsaturated profile is clarified, along with the key role of θ m .

IV. SATURATED ALIGNED AND MID-WAY FLUX LINKAGES
This section describes the idea proposed for the evaluation of the aligned and mid-way flux linkage curves: Ψ ph (i ph , θ a ) and Ψ ph (i ph , θ m ).Before embarking on the model description, it is highlighted that the formulation used in this work for the Ψ ph (i ph , θ) loci is based on the definition of a set of continuous functions Ψ ph vs. θ, each for a discrete value of i ph (as illustrated in Fig. 1).Hence, from now on, the following notation is applied: 1) Ψ ph ([i ph ], θ), denoting the flux linkage loci obtained for discrete values of i ph , 2) i * ph , a given value of i ph .According to this notation, the objective in this section is to determine the flux linkages for any i * ph with the rotor being "locked" respectively at θ a and θ m .The outlook of these curves can be seen in the red curves highlighted in Fig. 4.

A. Aligned Position
At full alignment, a typical flux density and flux lines distribution is shown in Fig. 7(a).Based on this distribution, it is widely accepted to model the magnetic circuit as shown in Fig. 8, where grey regions represent the portions of stator and rotor yokes having a non-negligible magnetic voltage drop, as yoke portions corresponding to unsupplied teeth offer very low reluctance due to the wide cross section [7].By assuming that the magnetic flux density is uniform within each of the above regions, Ampere and Gauss equations can be written as in ( 9)- (13), which have been expressed by taking the stator tooth cross section A st in (8) as a reference for the magnetic  flux crossing a magnetic pole, φ pole , neglecting the difference between tooth cross section and area of the tooth face facing the airgap for both stator and rotor (curvature effect).Besides, l sy and l ry represent the magnetically active path length of stator and rotor yoke (grey zones in Fig. 8), whose expressions can be found in [7].
The meaning of the geometrical quantities in ( 8)-( 13) is indicated in Table I At this point, ( 9)-( 13) can be solved for any given i * ph and the value of B st (i * ph ) is then easily obtained.By assuming that all coils of the same phase are series-connected, Ψ ph (i * ph , θ a ) may be obtained via (14), where P denotes the number of stator teeth pairs per phase:

B. Airgap Factor
The concept of airgap factor K g , relatively similar to the traditional Carter's factor, is introduced for taking into account the magnetic voltage drop caused by the fringing flux in the airgaps [31].By bearing in mind that the proposed model is mainly focused on early design stages, with the intention of maintaining a good compromise between model complexity and accuracy, the simplifying hypothesis is made that K g is not a function of the saturation level (an analysis is provided later in Section VII-B).Based on this hypothesis, K g (θ a ) can be evaluated in the following way.At the aligned position, Ψ ph can be evaluated through ( 8)-( 13), which, in absence of saturation, simplify into (15): On the other hand, at a low i * ph and hence no saturation, Ψ ph may be also obtained from L u a .Therefore, by plugging the value of B st obtained from (15) into ( 14), K g (θ a ) may be found by equalling Ψ ph to the product L u a • i * ph , resulting in: Another interpretation of K g (θ a ) can be given in relation to (13).Based on the Gauss' law, φ pole is a constant, so that for stator tooth and airgap, one finds: Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.where A eq (θ a ) denotes an equivalent airgap cross section that takes the fringing flux into account.Therefore, by observing (13), it is found that K g (θ a ) also represents the ratio between A st and A eq (θ a ), which can therefore be evaluated as follows: To conclude this subsection it is highlighted that the Carter factor is traditionally defined for singly-salient AC machines and is calculated from geometric parameters [31], [32].However, in the case of SRMs, double salience makes a pure geometric derivation far more complex.To this end, this work tackled this challenge by using the inductance, as in (16).

C. Mid-Way Position
An example of flux density and flux lines at mid-way position is shown in Fig. 7(b).As it can be noted, the main challenge in modelling this magnetic behaviour stems from the uneven flux density distribution within stator and rotor teeth, with flux lines that can no longer be approximated as straight lines.To provide a relatively simple approximation to this behaviour, this work proposes a discretisation of stator and rotor teeth in a suitably large number, n, of planar active slices conveying all the magnetic flux, as shown in Fig. 9.Further into detail, this simplification assumes all magnetic flux passing through the flux tube enclosed by the profile defined in (19), with a constant magnetic flux density inside each slice, whereas in real cases the magnetic flux crosses the entire cross section -including the teeth corner -and the magnetically equipotential surfaces are not planar.
In the case of the stator teeth, as one moves from the airgap towards the yoke, the i th slice possesses the following cross section A i and thickness l i : As it may be noted, the cross section is modelled via a hyperbolic tangent shape.This choice has been made due to the relative similarity with flux lines observed via FEA.In Fig. 7(b), a hyperbolic tangent profile is overlapped to the flux lines.As it can be seen, both profiles match relatively well, including in the lower tooth region.Moreover, the 1st segment of stator and rotor, along with the airgap, possess an equivalent cross section A eq (θ m ) that takes fringing flux effects into account.Its value is obtained through the airgap factor at the mid-way position K g (θ m ), whose definition is analogous to (16): For the rotor teeth, a similar discretisation is adopted by using expressions analogous to (20) and (19).Another point to note is related to the final segments corresponding to i = n in (19).Here, given that n is sufficiently large, the cross section is almost equal to A st .Besides, it is worth noting that the factor 2 in the hyperbolic tangent argument in (19) has been selected to optimise the accuracy of the proposed method, whilst minimising the computation time.At this point, Ampere's and Gauss' equations can be written as in ( 23)- (28), where A eq (θ m ) is taken as a reference for φ pole : Then, ( 23)-( 28) can be solved numerically for any given i * ph .In particular, B st( 1) in ( 28) denotes the magnetic flux density at the lowest stator tooth active slice.From its value, Ψ ph (i * ph , θ m ) may be obtained via (29):

V. NON-LINEAR FLUX LINKAGE LOCI MODEL
This section describes the idea proposed for the interpolation technique of the Ψ ph ([i ph ], θ m ) curves.Before starting with the model description, it shall be noticed that for the proposed interpolation process, rather than using directly the flux linkage values, it is preferred to pass through the incremental inductance values, as it simplifies the formulation.In particular, the idea is to define a set of inductance vs. θ curves for discrete currents L ph ([i ph ], θ m ), from which the evaluation of Ψ ph ([i ph ], θ m ) is straightforward.An example of saturated inductance curve is shown in Fig. 10 (yellow curve).

A. Saturated Inductance Profile: Region I
Even in saturated conditions, aiming to minimise the model complexity, it is commonly accepted to neglect the effects of fringing flux and rounding in this area.Hence, the inductance is constant and equal to the value at maximum alignment: It is observed that in (33), Ψ ph (i * ph , θ a ) represents the flux linkage at a given i * ph and θ a , evaluated in Section IV-A.

B. Saturated Inductance Profile: Region II
Region II represents the first half of the partial overlap region, i.e., between θ 1 and θ m .In this region, the magnetic behaviour is dominated by bulk saturation effects, producing relatively uniform saturation levels within teeth and yokes.This results in a nonlinear trend of the inductance curve at relatively high current levels, as represented in Fig. 10.In order to model this nonlinearity, the aforesaid 2nd-order FK ( 6) is set by imposing the boundary conditions in (31).

C. Saturated Inductance Profile: Region III
Region III represents the second half of the partial overlap region, spanning from θ m to θ 2 .In this region, the magnetic behaviour is dominated by local saturation effects, which arise mostly in the vicinity of the teeth corners.In real conditions, the passage from the bulk saturation typical of Region II to this new condition is progressive.However, in a bid to reflect this behaviour in a simplified mathematical way, the mid-way position θ m has been chosen to mark the passage from one behaviour to another.Then, the smoothness of the transition is ensured by imposing the continuity of the profile's derivative at either side of θ m .
The reason why θ 2 is used in-lieu of θ 2 is now explained.The absence of bulk saturation in Region III permits to approximate it through a linear trend (see Fig. 10).On the other hand, with relatively high currents, corner saturation effects are observable even at position θ 2 , as it can be seen in Fig. 10, where the saturated profile at position θ 2 stands below the unsaturated one (red dotted curve).To take this effect into account, Region III is prolonged up to a new rotor position, i.e., θ 2 (i * ph ).This last is found by imposing the following geometrical constraints, resulting in the closed-form expression (32) shown at the bottom of the next page: 1) The inductance profile of Region III must pass through the mid-way inductance point L(i * ph , θ m ), 2) Inductance profile is tangent to the unsaturated profile of Region IV (red dotted line), previously introduced in Section III.

D. Saturated Inductance Profile: Region IV
Once the position θ 2 (i * ph ) has been identified, the hypothesis of negligible saturation effects in the non-overlapping zone may be made.Hence, Region IV spans between θ 2 (i * ph ) and θ u following the unsaturated profile of Section III-C.
Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.Fig. 11.Comparison between the proposed loci and the one reported in [9] in the region between θ a and θ m .

E. Discussion About the Proposed Profile Formulation
This subsection is devoted to a brief comparison of the interpolation technique proposed in this work against that developed in [9], whose key points are as follows: 1) Only three regions are defined, as position θ 1 is ignored.
In other words, Region I proposed in this work does not exist, whereas only Region II is considered between θ a and θ m , 2) Both nonlinear regions (Region II and Region IV) are modelled by 1st-order FK equations.An example of profile proposed in this work is compared against that proposed in [9] in Fig. 11, together with FEA values.The example is taken from the SRM-A analysed in the next section at maximum rated current.Here, the first point to note is that a 2nd-order equation caters for a shape that is perfectly smooth and yet closer to reality.The second point to note is concerned with the use of four regions rather than three, as keeping Region I and Region II separated sensibly increases the accuracy nearby θ 1 .

A. SRM Designs Description
As already mentioned in the introduction, the proposed model has been developed with the objective of maintaining high accuracy for SRMs featuring different geometries and conceived for different rated operating conditions.To this end, four SRM designs have been considered, with the idea to cover different  applications of modern SRMs, namely automotive (SRM-A and SRM-D) and flywheel energy storage (SRM-B and SRM-C), and thus cover an extensive variety of geometric features, power and speed ratings and DC-bus voltages.From now on, designs will be referred to as SRM-A to SRM-D.The main geometric and control parameters are listed in Table II, whereas Fig. 12 provides an outlook of the cross sections.

B. Validation Process Description
Reminding that the analytical model has been conceived for early design stages, most of the comparisons and validation are conducted against FEA.Nonetheless, experimental results from SRM-A have been also included, with the purpose of validating the FEA results and assessing the discrepancy between Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

TABLE III FLUX LINKAGE ERRORS AT MID-WAY POSITION
analytical and experimental results.In its test bench shown in 13(b), SRM-A has been coupled to an induction machine through a Magtrol TM300 in-line torque transducer.Converter and control platform are also depicted, this last being realised on a DSP/FPGA platform [33].In the experimental testing, SRM-A has been run at a fixed speed and measuring current, voltage and average torque.Eventually, phase flux linkage and instantaneous torque have been obtained by post processing experimental voltages and currents (discrepancy between measured and post-processed average torque is 1.5%).More details about the FEA models, as well as about the testing bench can be found in [30].

C. Aligned and Mid-Way Flux Linkage Curves Validation
The flux linkage vs. current loci of the four designs at both aligned and mid-way rotor positions are shown in Fig. 14(a), (b), (c), and (d), where a very close similarity between analytical (solid lines) and FEA (dashed lines) results can be observed.Errors incurred at mid-way position are represented in Fig. 14(e) and are reported in Table III, where maximum, minimum and Mean Absolute Error (MAE) from amongst all loci points considered are shown (% errors are referred to the FEA values).As it can be seen, maximum and minimum errors are below 10%, with MAEs remaining below 5% for all designs.In addition, by observing Fig. 14(e), for all designs, the maximum positive error at θ m is incurred in correspondence of the knee-point of the saturation curve.This discrepancy stems from the hypothesis that all magnetic flux flows in a single tube, whereas in the real case this transition is smoothed by fringing and leakage phenomena.

D. Complete Flux Linkage Loci Validation
Flux linkage loci of the four designs are shown in Fig. 15

TABLE IV FLUX LINKAGE LOCI, MEAN ABSOLUTE ERRORS
plotted in the rotor position vs. current planes in Fig. 15(b),  (d), (f), and (h).In terms of error, Table IV reports the MAEs from amongst all loci points considered.As it can be noted, all values are lower than 4%.

E. SRM Performance Prediction Validation
In Table V, average torque, torque ripple and RMS currents obtained analytically, via FEA and experimentally for SRM-A are compared.As it is common practice for design processes, analytical and FEA performance have been obtained by considering a constant converter's DC-bus voltage and a hysteresis current control set to provide the desired torque.For the SRM-A, the same control parameters as in the physical prototype have been used.In terms of average torque, errors incurred are below 10% for all designs, even when analytical results are compared against experimental values.The flux linkage vs. current loops of the four designs are illustrated in Fig. 16(a), (c), (e), and (g).These last are reported as from their subtended area the Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.average torque per phase may be obtained, in accordance with the energy/co-energy theory (see (1)).
In terms of torque ripple, very good accuracy is obtained for SRM-A, SRM-B and SRM-C.For SRM-D, a 13.13% error is incurred.Instantaneous torque waveforms are shown in Fig. 16(b), (d), (f), and (h).Finally, regarding the RMS currents, errors incurred against FEA are below 4%.With respect to the experimental value, discrepancies around 10% are obtained for both analytical and FEA results.These are mostly caused by a constant voltage drop in the physical converter prototype not considered in analytical and FEA models, as it is demonstrated in [30].Finally, it is observed that the discrepancy between FEA and experimental flux linkage vs. current loops is caused by non-ideal aspects of physical prototypes, such as mutual coupling and uneven airgap, which create a "domino" effect when the hysteresis controller is engaged, as they cause switching instances to take place at different rotor positions in the FEA and experimental case, thereby causing the different trajectories.

A. Model Implementation in a Design Routine
This subsection provides a brief benchmark design exercise for SRM-A, [7], and SRM-B to SRM-D, with the objective of showing the reduction in computational time obtained by using the proposed model.To this end, a Genetic Algorithm is run, set on 50 generations with 15 design candidates each.The routine is implemented in MATLAB and runs on a workstation with an i7-3630 processor @2.40 GHz, 24 GB RAM.The computation times needed to obtain the Ψ ph ([i ph ], θ) curves of 750 candidates are compared in Table VI.As it can be seen, the fully analytical model allows one to complete the process in less than 1 h and is between 357 and 1785 times faster than the FEA process.

B. Saturation Effects on the Airgap Factor
Following from Section IV-B, the simplifying hypothesis of considering K g independent of saturation is now verified via FEA.For both θ a and θ m , in unsaturated and saturated conditions, a static simulation is run to evaluate the airgap flux density underneath one stator pole (spanning the entire angle ζ s ), from which φ pole−F EA is found.Subsequently, in accordance with the literature, [31], [32], the peak flux density B g−F EA is used to attain A eq−F EA in accordance with (17) (B g in (13) needs to be interpreted here as a peak value): Finally, K g−F EA is derived via (18), where A st and A st /2 are considered for θ a and θ m respectively.Results are shown in Table VII, where a considerably small variation is observed between unsaturated and saturated conditions, confirming the viability of the bespoke hypothesis.

C. Results Comparison With First-Order FK Equation
In order to provide an example of the improvement introduced by the proposed model against that proposed in [9], this last is implemented for SRM-D to predict the flux linkage and torque characteristics.Results are shown in Fig. 17 and are quantified in Table VIII, where the improvement attained with the proposed model can be noted.

D. Model Implementation in High-Speed Conditions
When designing for high-speed conditions, different considerations should be made.In fact, SRMs operate in the socalled 'single-pulse' mode, with the highest phase currents being reached in poor overlap conditions, so that only local saturation in the teeth corners arises.As a result, even analytical models Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
that neglect the effects of saturation at high speed may ensure acceptable accuracy for early design stages [21].

E. Model Implementation With Tapered-Shaped Teeth
In modern SRMs, tapered teeth where side edges are not straight are not unusual [33].In this case, the assumption of modelling the flux tube at θ m via a hyperbolic tangent needs to be verified, as higher estimation errors might be incurred.

VIII. CONCLUSION
This work proposed a novel design-oriented analytical model predicting the Ψ ph (i ph , θ) loci for SRMs of any geometry and rated operating conditions.The model relies on a new interpolation technique based on 2nd-order FK equations, as well as on an analytical model catering for the flux linkage in partial overlap and saturated conditions.Further to the validation against four FEA models and one physical prototype, the following key conclusions can be drawn: r Use of 2nd-order FK equations in lieu of 1st-order ones leads to a perfectly smooth yet more accurate representation of Ψ ph curve, with considerable improvements in SRM performance prediction.
r Very high accuracy in predicting Ψ ph in partial overlap and saturated conditions is obtained, with MAEs lower than 4%.
r Very high accuracy in predicting the main SRM perfor- mance is also obtained, with errors lower than 10% against FEA.r A reduction in computation time of between 357 and 1785 times compared to the FEA process is attained.
r An alternative yet useful formulation of the airgap factor has been proposed based on the inductance value, whose independence of saturation is demonstrated.

Fig. 3 .
Fig. 3. High-Level overview of the analytical SRM performance prediction process used in this work.

Fig. 4 .
Fig. 4. Schematic example of flux linkage vs. phase current loci, plotted for constant rotor positions.

Fig. 5 .
Fig. 5. Example of inductance vs rotor position curve in unsaturated conditions.

Fig. 7 .
Fig. 7. Magnetic Flux Density and Flux Lines distribution at the aligned position (a) and mid-way position (b).
, whereas new quantities introduced are as follows: N : Number of series connected turns per stator teeth pair, B st : Stator tooth magnetic flux density, B rt : Rotor tooth magnetic flux density, B g : Mid-Airgap magnetic flux density, B sy : Stator yoke magnetic flux density, B ry : Rotor yoke magnetic flux density, K g (θ a ): Airgap factor at aligned position, μ i (B i ): Magnetic permeability of the i th circuit part, l sy : Stator yoke active path length, l ry : Rotor yoke active path length.

Fig. 10 .
Fig. 10.Example of inductance vs rotor position curve in unsaturated and saturated conditions.
31) In particular, the fourth boundary condition sets the derivative equal to the slope of the profile in Region III.Here, two observations are needed: 1) Rotor position θ 2 (i * ph ) is different from θ 2 and is a function of i * ph , as discussed in the next subsection, 2) Ψ ph (i * ph , θ m ) appearing in the second equation of (31) represents the flux linkage at a given i * ph and θ m , evaluated in Section IV-B.

Fig. 12 .
Fig. 12. Geometry outlook of the SRM designs used for model validation.

TABLE I GEOMETRIC
PARAMETERS OF AN SRM DESIGN CANDIDATE

TABLE II GEOMETRIC
PARAMETERS OF THE SRM DESIGNS

TABLE VIII FIRST
-ORDER FK EQUATION RESULTS FOR SRM-D