Multi-Physics Multi-Objective Optimal Design of Bearingless Switched Reluctance Motor Based on Finite-Element Method

Abstract

The bearingless switched reluctance motor (BSRM) integrates the switched reluctance motor (SRM) with the magnetic bearings, which avoids mechanical bearings-loss and makes it promising in high-speed applications. In this paper, a comprehensive framework for the multi-physics multi-objective optimal design of BSRMs based on finite-element method (FEM) is proposed. At first, the 2-D electromagnetic model of a fabricated initial design prototype is built and solved by the open-source FEM software, Elmer. The iron loss model in Elmer based on the Fourier series is modified by a transient iron loss model with less computation time. Besides, a simplified lumped-parameter (LP) thermal model of the BSRM is applied to estimate the temperature rise of BSRM in the steady state. Then, the comprehensive framework for the multi-physics multi-objective optimal design of BSRMs based on FEM is proposed. The objectives, constraints, and decision variables for optimization are determined. The multi-objective genetic particle swarm optimizer is utilized to obtain the Pareto front of optimization. The electromagnetic performance of the final optimal design is compared with the initial design. Comparison results show that the average electromagnetic torque and the efficiency are significantly enhanced.

1. Introduction

Bearingless switched reluctance motors (BSRMs) have the same structure as iron core with the conventional switched reluctance motors (SRMs), hence inheriting their advantages, e.g., low cost, robust rotor structure, and fault tolerance capability. In addition, the BSRM integrates the SRM with the magnetic bearings, which avoids mechanical bearings loss and makes it promising in high-speed applications, i.e., flywheel energy storage [1,2,3]. Compared with the SRMs, the rotation and suspension of the rotor in BSRMs should be controlled simultaneously. Control strategies for the BSRMs with dual-winding and single windings have been investigated by many researchers [2,3,4,5,6,7,8,9,10]. Besides, the BSRMs with novel topologies were also presented [11,12,13]. However, the studies on the optimal design of the BSRMs are limited. At present, the optimal arrangement of windings was discussed in detail [14,15]. The multi-objective optimization of the dual-winding BSRMs based on the analytical model was presented [16]. Nevertheless, the analytical electromagnetic model is simplified and the magnetic saturation is unconsidered. The detail dimensions of stator and rotor poles are not optimized, and the temperature rise of the BSRM cannot be considered during the optimization.

In order to consider the magnetic saturation, the non-linear models based on artificial neural networks (ANNs) or response surface (RS) have been proposed [17,18]. However, the design of the ANNs or RS model is complicated. The built ANNs or RS model also requires additional validation of its accuracy in the objective space of optimization. Compared with the analytical models, the non-linear electromagnetic characteristics of electrical machines with complicated shape can be solved accurately by the finite-element method (FEM), which makes it more popular for the analysis and optimal design of the electrical machines [19,20,21]. Nevertheless, the commercial FEM software are always integrated and encapsulated only with application programming interface (API) [22]. The detail algorithms of electromagnetic solver used in the commercial software are unknown and difficult for users to change or modify the functions. Hence, some open-source softwares become remarkable. For instance, the open-source FEM software, Elmer, has been used to model the skewed electrical machine via a multi-slice model [23].

In addition, the optimal design based on the multi-physics models of electrical machines has become a research hotspot [24,25,26,27]. The reason is that the temperature rise in the windings determines the insulation level and service life of the electrical machine. Besides, the computational fluid dynamics (CFD) analysis is essential for the establishment of an accurate thermal model for the prediction of the temperature rise. However, more time consumptions are required for the CFD analysis.

In this paper, a multi-physics multi-objective optimal design of the BSRM based on FEM is proposed. At first, the principle of the dual-winding BSRM is briefly introduced. The electromagnetic model is established by the open-source FEM software package, Elmer [23]. In order to decrease the time consumption for calculating the iron loss of laminations, a transient iron loss model modifies the iron loss model in Elmer. Estimation results of iron loss are analyzed and compared. Second, a simplified lumped-parameter (LP) thermal model of the BSRM is used to estimate the thermal performance of BSRM in the steady state. Due to the fact that the surface velocity of the rotor is not high enough, the structural analysis of rotor is unconsidered. Third, the objective functions, constraints, and decision variables of optimization of BSRMs are determined. The multi-physics model of BSRMs based on the Elmer and LP thermal model is integrated and the comprehensive framework for multi-physics multi-objective optimal design of BSRMs is proposed. Corresponding parallel computation for optimization is also investigated. At last, the optimization results, i.e., the solutions of the searched Pareto front, are obtained and analyzed. The final optimal design is chosen from the Pareto front. Its performance is compared with that of the initial design. Comparison results illustrate the effectiveness of the proposed optimal design method.

2. Finite-Element Model of BSRMs

2.1. Principle of Dual-Winding BSRM

The winding configuration of the A-phase for the three-phase 12/8-pole BSRM is shown in Figure 1. W_ma represents the main windings on the stator poles. The suspension windings in α- and β-direction are denoted as W_sa1 and W_sa2, respectively. When the current of main windings i_ma and current of suspension windings i_sa1 are excited with the same direction shown in Figure 1, the flux density in the air-gap a₁ would superimposes and the flux density in the air-gap a₃ would be weaken. Then, a radial suspension force towards the α-direction can be generated. Similar principle for generating the suspension force at the β-direction can also be obtained. Hence, the position of rotor in two perpendicular directions can be controlled by regulating the current of main and suspension windings [2].

2.2. 2-D Finite-Element Model of BSRM

The three-phase 12/8-pole BSRMs use dual-windings, i.e., the main and suspension windings, to generated different magnetic flux density in the air-gap. Then, the radial force is produced and acts on the rotor for suspension [2,3]. The mathematical models of suspension force and electromagnetic torque have been investigated [2,5]. However, the model is linear and the permeance in the air-gap is simplified, which neglects the influences of detail dimensions, e.g., the width of stator or rotor poles, etc. Thus, in this paper, the 2-D FEM model is used to optimize these detail dimensions.

Elmer is an open-source FEM software package, which can solve partial differential equations for multi-physics simulations [23,28]. In Elmer, based on the Maxwell equations, the 2-D magnetic system is described as [29,30]:

$$\sigma \left(\frac{\partial }{\partial t}{A}_Z+\nabla V\right)+\nabla \times \left(\frac{1}{\mu }\nabla \times A_Z\right)=J_Z$$

(1)

where σ is the electrical conductivity, A_Z is the out-of-plane component of vector potential, V is the electrical scalar potential which is assumed to zero because the model is input with the current source, μ is the permeability of material, J_Z is the current density in the out-of-plane direction. After (1) is solved with proper boundary and initial conditions, the flux density B can be obtained by solving:

$$B=\nabla \times A_Z$$

(2)

Before the optimization of the dual-winding BSRMs, an initial design prototype was fabricated. The power and speed of the prototype are 1.5 kW and 5000 rpm, respectively. The DC-bus voltage of the power converter is 150 V. The sketch of the cross-section of laminations and corresponding parameters of the initial design are shown in Figure 2 and

Table 1. Main parameters of the initial design prototype.

Parameters	Value	Parameters	Value
r	25 mm	δ	0.3 mm
Nm	13	Ns	25
bs	6.6 mm	br	6.6 mm
hs	29.3 mm	hr	7 mm

, respectively. r is the radius of rotor pole, δ is the length of air-gap. N_m and N_s are the number of turns of the main windings and suspension windings, respectively. b_s and b_r are the width of stator pole and the width of rotor pole, respectively. h_s and h_r are the height of stator pole and the height of rotor pole, respectively. h_sy and h_ry are the yoke thickness of stator and the yoke thickness of rotor, respectively. d_shaft is the diameter of the shaft. The 0.5 mm silicon steel 50 W350 is used as the material of iron core. The magnetization curve is shown in Figure 3. The 2-D geometry model of the machine and corresponding mesh are built in Gmsh [31]. The mesh in the air-gap and top of the stator and rotor poles are dense enough, as shown in Figure 4. Then, the ElmerGrid, which is a part of Elmer package, is used to transform the mesh files and unite them to one mesh model [32]. At last, the material data, excitations, boundary conditions, and solvers are configured in the solver command file of Elmer.

If the electromagnetic model is solved by the time-stepping FEM with 1° mechanical angle per time-step, for exciting three phases in turn and 15° for each phase, 46 time steps are required. In order to decrease the computational effort during the optimization, the time-step is set as 3° and 16 time steps are calculated. The obtained electromagnetic torque and suspension force of the vertical direction using different time-step are shown in Figure 5a,b, respectively. For the time-stepping FEM with 3° time-step, in each phase, the values at four angular positions are used to estimate the average torque and the average suspension force. The calculated average torque with the 3° time-step and that with 1° time-step are 0.1796 Nm and 0.1623 Nm, respectively. Besides, the calculated average suspension force with the 3° time-step and that with 1° time-step are 45.98 N and 45.85 N, respectively. Corresponding percentage error is approximate 0.2835%. Thus, the results obtained with 3° time-step are almost identical as those with 1° time-step. The solved flux density can be visualized by the open-source post-process software ParaView [33], as shown in Figure 6. The values of superposed and counteracted flux density in the air-gap are approximate 0.6372 T and 0.01528 T, respectively.

2.3. Transient Iron Loss Model

During the optimization of the BSRMs, the iron loss of laminations should be considered for calculating the electrical machine efficiency. In Elmer, the Steinmetz equations with Fourier series are used to estimate the iron loss [30]. However, the distribution of flux density in different parts of laminations for SRMs or BSRMs is complicated. The calculation of the average iron loss based on the Fourier series requires at least one round of the rotor, which needs more computation time during the optimization. Besides, for SRMs or BSRMs, the high frequency harmonic components of the flux density exist in the stator or rotor core. However, only limited harmonics can be considered by the Fourier series.

Hence, the source code of the solver in Elmer is modified by a transient iron loss [34]. The iron loss per unit volume P_{iron_V} is estimated by [34]:

$$\begin{array}{ll}{P}_{\mathrm{iron}\_\mathrm{V}}& =P_{\mathrm{hys}\_\mathrm{V}}+P_{\mathrm{eddy}\_\mathrm{V}}\\ & =\left(H_{\mathrm{c}}+k_{\mathrm{h}}\left|B\left(t\right)\right|\right)\left|\frac{dB\left(t\right)}{dt}\right|+k_{\mathrm{eddy}}{\left(\frac{dB\left(t\right)}{dt}\right)}^2\end{array}$$

(3)

where P_{hys_V} and P_{eddy_V} denote the hysteresis loss and the eddy loss per unit volume, respectively. H_c, k_h, and k_eddy are model coefficient and can be fitted by the iron loss data of lamination provided by the manufacturer. The flux density B and dB/dt of each element can be obtained based on the solved electromagnetic field of the time-stepping FEM. In this paper, 16 time steps with the 3° time-step size of the time-stepping FEM are used to estimate the iron loss of the initial design prototype under 5000 rpm. The calculated average iron loss per unit volume based on the Fourier series is shown in Figure 7a. The iron loss in the stator is non-uniformly distributed. The iron loss is also greater at the top corners of the rotor poles. This phenomenon is caused by the greater dB/dt at the corners. Besides, the transient distribution of the iron loss at 0° calculated by the transient iron loss model is shown in Figure 7b, which basically agrees with results shown in Figure 7a. In addition, the estimated iron loss of the initial design prototype based on the Fourier series and the transient iron loss model are listed in

Table 2. Comparison of the iron loss in different parts.

Parts	Fourier Series Considering 25 Harmonics (W)	Fourier Series Considering 36 Harmonics (W)	Transient Iron Loss Model (W)
Rotor	3.926	4.872	6.255
Stator	19.32	23.82	27.71
Total iron loss	23.24	28.69	33.97

. The first 25 and 36 harmonic components of flux density are considered by the Fourier series, respectively. Table 2 reveals that considering more harmonics can lead the iron loss calculated by the Fourier series approaching the results of the transient model. However, more calculations are required if more harmonic components are considered. Hence, in this paper, the transient model is used to estimate the iron loss for optimizations.

3. Thermal Model of BSRMs Based on Simplified LP Thermal Networks

The maximum temperature rise of the electrical machine directly influences the insulation level and the service life. Thus, the estimation of the temperatures in different regions in the machine is necessary. In this paper, the temperature in the steady-state is used to evaluate the thermal performance of the BSRMs. Thus, the simplified LP thermal model of BSRMs is applied, as shown in Figure 8. The thermal capacitances are neglected. The basic elements and other equivalent thermal resistances form the LP thermal network. Γ_amb is the temperature of ambient. P_sy, P_sp, P_rp, P_ry and P_copper are the iron loss of the stator yoke, the iron loss of the stator pole, the iron loss of the rotor pole, the iron loss of the rotor yoke and the copper loss of windings, respectively. R_{sy_1r}, R_{sy_2r} and R_{sy_3r} are the radial thermal resistance of stator yoke. R_{sp_1r}, R_{sp_2r} and R_{sp_3r} are the radial thermal resistance of stator pole. R_{rp_1r}, R_{rp_2r} and R_{rp_3r} are the radial thermal resistance of rotor pole. R_{ry_1r}, R_{ry_2r} and R_{ry_3r} are the radial thermal resistance of rotor yoke. The hollow cylindrical element is used to describe the heat transfer in the poles and yoke of stator and rotor [35,36,37]. Corresponding thermal resistances can be calculated by the geometrical dimensions of the lamination [37]. R_{case_amb} is the natural convection thermal resistance between case and ambient. R_{cap_ia} is the convection thermal resistance between endcaps and inner air. R_{case_win} is the convection thermal resistance between case and end windings. R_{ia_win} is the convection thermal resistance between inner air and windings. R_insul is the conduction thermal resistance of insulation. R_{ry_shaft} is the conduction thermal resistance between rotor yoke and shaft. R_{sp_air} and R_{air_rp} are the convection thermal resistance between stator inner surface and air-gap and the convection thermal resistance between air-gap and rotor outer surface, respectively. R_{sy_case} is the conduction thermal resistance between stator yoke and case. For these thermal resistances, the methods to determine the value of R_{cap_ia}, R_{case_win}, R_{ia_win}, R_{ry_shaft} have been presented [38]. In this paper, for the salient pole rotor of BSRM, the CFD is utilized to determine the value of R_{sp_air} and R_{air_rp}. The thermal resistances, i.e., R_{case_amb}, R_{sy_case}, R_insul, are calibrated by direct current thermal tests [38].

4. Comprehensive Framework for Multi-Physics Multi-Objective Optimization of BSRMs

4.1. Objective Functions for Optimal Design of BSRM

As the electrical machine, maximizing the electromagnetic torque of the BSRM with the least cost of iron and copper materials is always the demand for the optimal design. Besides, the efficiency is another important criterion. In this paper, two objectives for optimizations are determined. The first objective is the specific torque T_m:

$$T_{\mathrm{m}}=\frac{T_{\mathrm{avg}}}{m_{\mathrm{iron}}+m_{\mathrm{copper}}}$$

(4)

where T_avg is the average torque of BSRM. m_iron and m_copper are the mass of the used iron and copper. The second objective is the efficiency of the BSRM η:

$$\eta =\frac{\omega T_{\mathrm{avg}}}{\omega T_{\mathrm{avg}}+P_{\mathrm{iron}}+P_{\mathrm{copper}}}$$

(5)

where ω is the angular velocity of rotor, P_iron and P_copper are the iron loss and copper loss, respectively. The average electromagnetic torque and iron loss of BSRM can be obtained by the time-stepping method in FEM model. It should be pointed out that these objectives have mutual influences and a tradeoff solution is required. Thus, the multi-objective optimization algorithm should be applied to obtain the Pareto front of the optimal design.

4.2. Constraints for Optimal Design of BSRM

Different from the conventional SRMs, the design of BSRM should consider the suspension force to guarantee the suspension of the rotor and shaft. Due to the fact that the suspension force changes with the angular position of the rotor, the average suspension force is used to evaluate the suspension force of the BSRMs. Considering the margin of suspension force, following constraint should be satisfied:

$$F_{\mathrm{avg}}\ge 2\left(m_{\mathrm{rotor}}+m_{\mathrm{shaft}}\right)g$$

(6)

where F_avg is the average suspension force, m_rotor and m_shaft are the mass of rotor and shaft, respectively. g is the acceleration of gravity.

In order to guarantee the safe operation of the BSRM, the maximum temperature in the steady-state should be limited. In this paper, considering the insulation level of the electrical machine, the maximum temperature of the BSRM should be less than 130 °C. However, considering the estimation error of the thermal model and overload performance of BSRM, the maximum temperature of the optimized BSRM Γ_max should satisfy:

$${\Gamma }_{\max }\le 0.8\times 130°\mathrm{C}$$

(7)

In addition, the filling factor of windings in the stator slots k_u should be constrained [16]:

$$k_{\mathrm{u}}\le 0.45$$

(8)

4.3. Decision Variables for Optimal Design of BSRM

In this paper, the outer diameter of stator core, the stack length, and length of air gap are kept as same as that of the initial design prototype. The power of the main windings and suspension windings should be distributed reasonably. The ratio between the power of the main windings and that of the suspension windings k_p is described as:

$$\{\begin{array}{l}{k}_{\mathrm{p}}=\frac{P_{\mathrm{m}}}{P_{\mathrm{s}}}\\ P=P_{\mathrm{m}}+P_{\mathrm{s}}\end{array}$$

(9)

where P_m and P_s are the power of main and suspension windings, respectively. P is the rated power of the BSRM. Then, the rated current of main and suspension windings, i.e., i_m and i_s, can be determined by the power of the main and suspension windings with the bus voltage U [16]:

$$\{\begin{array}{l}{i}_{\mathrm{m}}=\frac{k_{\mathrm{I}}{i}_{\mathrm{mp}}}{\sqrt{3}}=\frac{k_{\mathrm{I}}}{\sqrt{3}}\frac{P_{\mathrm{m}}\left({\eta }_{\mathrm{a}}+1\right)}{2{\eta }_{\mathrm{a}}U}\\ i_{\mathrm{s}}=\frac{k_{\mathrm{I}}{i}_{\mathrm{sp}}}{\sqrt{3}}=\frac{k_{\mathrm{I}}}{\sqrt{3}}\frac{P_{\mathrm{s}}\left({\eta }_{\mathrm{a}}+1\right)}{2{\eta }_{\mathrm{a}}U}\end{array}$$

(10)

where k_I is used to consider the margin of the current in each phase, η_a is the assumed efficiency [16]. Besides, the total flux density in the air-gap is superposed by the flux densities independently generated by the main and suspension windings. The total flux density in the air-gap where the flux density is superposed B_total is expressed as:

$$B_{\mathrm{total}}=B_{\mathrm{m}}+B_{\mathrm{s}}$$

(11)

where B_m and B_s are the flux density independently generated by the main and suspension windings, respectively. Thus, the ratio between the flux density generated by the main windings and that generated by the suspension windings k_b is described as:

$$k_{\mathrm{b}}=\frac{B_{\mathrm{m}}}{B_{\mathrm{s}}}$$

(12)

The power ratio k_p, ratio of flux density k_b, and the total flux density B_total should be optimized. Then, the number of turns of main windings and suspension windings can be estimated as [16]:

$$\{\begin{array}{l}{N}_{\mathrm{m}}{i}_{\mathrm{m}}\approx \frac{B_{\mathrm{m}}{k}_{\mathrm{m}}\delta }{{\mu }_0}\\ N_{\mathrm{s}}{i}_{\mathrm{s}}\approx \frac{B_{\mathrm{s}}{k}_{\mathrm{m}}\delta }{{\mu }_0}\end{array}$$

(13)

where k_m is a coefficient for compensating the magnetic potential difference in the iron core [16]. The geometrical parameters, the inner diameter of the stator D, the pole arc coefficient of stator α_s, and the pole arc coefficient of rotor α_r, are optimized. The length of air-gap, the thickness of stator, and rotor yoke are kept as constants. All the decision variables and constants are listed in

Table 3. Decision variables and constraints for optimization.

Parameters	Type	Value
kp	Variable	2–10
Btotal	Variable	0.4–1.4 T
kb	Variable	2–20
D	Variable	40–65 mm
αs	Variable	14°/30°–20°/30°
αr	Variable	14°/45°–20°/45°
δ	Constants	0.3 mm
hsy	Constants	5.4 mm
hry	Constants	7 mm
dshaft	Constants	22 mm
L	Constants	80 mm

. The pole arc coefficient is defined as the ratio between the pole arc and pole pitch.

4.4. Comprehensive Framework for Multi-Physics Multi-Objective Optimization of BSRMs

In this paper, a multi-objective optimizer is necessary. The multi-objective genetic particle swarm optimizer (MOGPSO) has been presented to solve the nonconvex problems with non-uniformly distributed Pareto front [16]. Thus, in this paper, the MOGPSO is utilized to search the Pareto front.

The proposed comprehensive framework for multi-physics multi-objective optimal design of the BSRMs is shown in Figure 9. The MOGPSO and the simplified LP thermal model are built in MATLAB. The position of each particle of the MOGPSO represents a candidate solution of optimal design. After initializing the parameters of the optimizer, the calculated dimensions of the stator and rotor are used to build the geometrical file and input to the Gmsh to obtain the original mesh files. Then, the ElmerGird is used to generate the mesh files for ElmerSolver. The definitions of excitations, boundary conditions, etc., are configured in the solver command file. The electromagnetic torque, suspension force and the iron loss are calculated by the time-stepping FEM. Then, the iron loss and copper loss data are used as the heat sources for the simplified LP thermal model to obtain the maximum temperature rise and its location in the BSRM. Thus, the objective functions of each particle can be calculated, and corresponding constraints are evaluated. In the MOGPSO, the concept of Pareto domination is applied to select the particles with better performance according to the values of two objective functions. The domination relationships between particles are determined and all the non-dominated particles are stored in the repository. Then, the crossover and mutation operators are used and generate crossover and mutation offspring. After determining the domination relationship between the particles in repository, crossover offspring, and mutation offspring, the non-dominated particles are filtered and stored in the repository. The group leader is selected based on the divided grid and roulette wheel selection. Then, the next iteration begins and the positions of all particles are updated as [16]:

$$\begin{array}{ll}Ve{l}_i\left(k+1\right)& =wVe{l}_i\left(k\right)\\ & +c_1{r}_1\left(k\right)\left(Po{s}_{i,\mathrm{pbest}}\left(k\right)-Po{s}_i\left(k\right)\right)\\ & +c_2{r}_2\left(k\right)\left(Po{s}_{\mathrm{leader}}\left(k\right)-Po{s}_i\left(k\right)\right)\end{array}$$

(14)

$$Po{s}_i\left(k+1\right)=Po{s}_i\left(k\right)+Ve{l}_i\left(k\right)$$

(15)

where Vel_i(k) and Pos_i(k) denote the velocity and position of the ith particle at the kth iteration; w, c₁, and c₂ represent the inertia weight and the learning coefficient of the individual best and the group leader, respectively; r₁(k) and r₂(k) are the uniformly distributed random numbers in the range of 0 to 1 at the kth iteration; Pos_i,pbest(k) and Pos_leader(k) are the position of the individual best of the ith particle and the group leader, respectively. After iterations, the positions of the particles would move towards the superior nondominated solutions and the Pareto front would converge. Corresponding flow chart of the MOGPSO for the multi-physics multi-objective optimal design of BSRMs is shown in Figure 10.

4.5. Parallel Computation of the Multi-Physics Multi-Objective Optimal Design of BSRMs

For the proposed multi-physics multi-objective optimization of BSRMs, the parallel computation is necessary to investigate. MATLAB provides the parallel for-loops parfor, which can be executed in parallel [39]. Hence, in this paper for the MOGPSO, the candidate solutions expressed by different particles are modeled and solved by different CPU cores. Then, more than one time-stepping FEM model can be solved simultaneously. In this paper, a DELL Precision Tower Workstation (T5810) is utilized for parallel computation of the proposed multi-physics multi-objective optimal design of BSRMs.

5. Optimization Results of Multi-Physics Multi-Objective Optimization of BSRMs

5.1. Analysis of Searched Pareto Front

In this paper, 48 particles are iterated 180 times to obtain the Pareto front of the multi-physics multi-objective optimal design of BSRMs. The searched Pareto fronts with the interval of 25 iterations are shown in Figure 11. At the beginning of the iterations, few solutions are searched. The reason is that the iteration is not convergent and the searched Pareto front is dominated by the searched Pareto front in the next iteration. After 175 iterations, the Pareto front basically converges. The finally searched solutions distributes evenly on the Pareto front for the 180th iteration.

Figure 12 reveals the convergence trend of the main parameters during the optimization. The average electromagnetic torque gradually increases. The calculated average suspension forces converge in the range of 46.42 N to 50.83 N and satisfy the constraint. Nevertheless, the stator inner diameters of most solutions in the final Pareto front converge to the same value. Similar results are obtained for the width of rotor poles, which are directly influenced by the decision variables α_r. The width of stator poles converge in the range of 6.357 mm to 7.092 mm. This result indicates that the dimensions of the laminations are basically convergent for the same solutions in the final Pareto front. The difference of the solutions in the final Pareto front is mainly caused by different number of turns of the suspension windings. It causes different values of suspension force, copper loss, and the maximum temperature of the BSRMs. The maximum temperatures of solutions in the Pareto front approach to the limitation of the constraint, i.e., 104 °C.

5.2. Electromagnetic and Thermal Performance of Final Optimal Design

After the Pareto front is obtained, the final optimal design should be determined to obtain the tradeoff solution for two objectives. A common method is to find the solution which has the minimum distance to the ideal solution in the objective space, as shown in Figure 13. The parameters of the lamination for the final optimal design are shown in

Table 4. Main parameters of the final optimal design.

Parameters	Value	Parameters	Value
r	22.38 mm	δ	0.3 mm
Nm	14	Ns	63
bs	6.636 mm	br	6.314 mm
hs	31.92 mm	hr	4.383 mm
im	8.250 A	is	0.8436 A

. The radius of the rotor poles is less than that of the initial design. Besides, the width of the stator poles is greater than the width of the rotor poles. This approach can increase the effective angular range of the electromagnetic torque and suspension force.

The calculated electromagnetic torques in A-phase T_a of the final optimal design and that of the initial design are compared, as shown in Figure 14a. The electromagnetic torque is enhanced. The generated suspension forces in the vertical direction F_β are shown in Figure 14b. The suspension force is slightly greater than that of the initial design. Figure 15 reveals the distribution of the flux density for the final optimal design. Compared with Figure 6, the superposed and counteracted flux density in the air-gap both increase accordingly. The values of superposed and counteracted flux density in the air-gap are approximate 0.7037 T and 0.2589 T, respectively. The detail comparisons of the performance are listed in

Table 5. Comparisons of enhancements of machine performance.

Parameters	Final Optimal Design	Initial Design	Compared Percentage Ratio (%)
Tavg	0.2045 Nm	0.1623 Nm	+26.00
miron	3.253 kg	3.328 kg	−2.254
mcopper	0.8234 kg	0.4428 kg	+85.95
Tm	0.0502 Nm/kg	0.0430 Nm/kg	+16.74
Piron	36.28 W	33.97 W	+6.800
Pcopper	9.181 W	6.177 W	+48.63
η	0.7020	0.6792	+3.357
Favg	46.77 N	45.85 N	+2.007
Γmax	103.97 °C	95.81 °C	+8.517
Location of Maximum Temperature	Windings	Windings	\

. The mass of the utilized iron has a minor difference between the final optimal design and the initial design. However, the average electromagnetic torque is enhanced. This enhancement is attributed to the reasonable distribution of the power and generated flux density between the main and suspension windings, which directly influences the rated current and the number of turns for the main and suspension windings. Corresponding iron loss and copper loss also slightly increase. However, the final optimal design obtains better performance and efficiency.

6. Conclusions

In this paper, the multi-physics multi-objective optimal design of the BSRMs based on FEM is proposed. First, the electromagnetic FEM model is built and solved. The iron loss model is modified by the transient iron loss model. Then, the simplified LP thermal model of the BSRM is used to estimate the maximum temperature rise in the steady-state. Then, the objective functions, constraints, and decision variables are determined and the comprehensive framework of the multi-physics multi-objective optimal design of BSRMs is proposed. The final optimal design is selected from the obtained Pareto front. The electromagnetic torque of the final optimal design is enhanced significantly and corresponding efficiency is increased. Further research would focus on the optimization of the BSRMs, which integrates control strategies to alleviate the torque ripple.