Analysis of a Shaftless Semi-Hard Magnetic Material Flywheel on Radial Hysteresis Self-Bearing Drives

Abstract

Flywheel Energy Storage Systems are interesting solutions for energy storage, featuring advantageous characteristics when compared to other technologies. This has motivated research effort focusing mainly on cost aspects, system reliability and energy density improvement. In this context, a novel shaftless outer-rotor layout is proposed. It features a semi-hard magnetic FeCrCo 48/5 rotor coupled with two bearingless hysteresis drives. The novelty lies in the use of the semi-hard magnetic material, lending the proposed layout advantageous features thanks to its elevated mechanical strength and magnetic properties that enable the use of bearingless hysteresis drives. The paper presents a study of the proposed layout and an assessment of its energetic features. It also focuses on the modeling of the radial magnetic suspension, where the electromagnets providing the levitating forces are modeled through a one-dimensional approach. The Jiles–Atherton model is used to describe the magnetic hysteresis of the rotor material. The proposed flywheel features a mass of 61.2 kg, a storage capability of 600 Wh at the maximum speed of 18,000 rpm and achieves an energy density of 9.8 Wh/kg. The performance of the magnetic suspension is demonstrated to be satisfactory and the influence of the hysteresis of the rotor material is highlighted.

1. Introduction

Nowadays, energy storage plays an essential role in the modern energy supply chain, in public transportation, urban grids and transport industry. It enables voltage sag compensation, frequency regulation, load peak shaving and power smoothing, thus improving the grid stability and optimizing the energy supply process. Additionally, it is fundamental to the mainstream adoption of renewable energy resources [1,2]. Flywheel Energy Storage Systems (FESSs) are devices that convert electricity to kinetic energy stored in a rotating mass. The flywheel is coupled with an electric machine acting as a motor or a generator during the charging or the discharging phase, respectively. A bearing system is needed to support the spinning flywheel. Active magnetic levitation is usually employed in modern systems since it features bearing losses lower than in conventional ball or journal bearings, lubrication-free operation, low maintenance costs along with the possibility to control the dynamics of the levitating body [3,4,5]. Research on active magnetic bearings (AMBs) has been focused on control and hardware improvements [6,7,8,9,10], system identification [11,12,13] and on integration between bearing and electric machine operation in the so-called self-bearing drive [14,15,16]. This configuration features lower component count and more compact layout than common solutions with separate electric motor and bearing system.

FESSs provide various advantages compared to other energy storage technologies, including high power rating, long lifetime, fast recharging, short response time and high power density. These devices require little periodic maintenance, have a very low environmental impact, their charge/discharge cycles do not degrade storage capacity, and their state of charge is easily obtained by measuring the rotational speed [17]. The use of self-bearing drives in flywheel applications lends the system more compactness, higher reliability, lower component count and mass saving, thus strengthening the FESS advantages [18,19,20,21]. On the other hand, one of the main drawbacks of FESSs is the low energy density when compared to other storage technologies [1]. This aspect has motivated a research effort devoted to investigating higher energy density layouts. One solution is represented by flywheels made of composite materials, which achieve higher energy density compared to high-strength steel flywheels thanks to their light-weight, along with the very large mechanical strength permitting very high rotational speeds. A rim made of composite material can be much lighter than a steel flywheel when storing a comparable amount of energy. In this configuration, the rotor is usually in the shape of a hollow cylinder made of composite material such as carbon fiber reinforced plastic.

Since the energy density depends directly on the mass of the whole FESS, it can be improved by designing as much mass of the system to contribute to storing energy. The index allowing quantification of this contribution is known as the Degree of Integration (DoI). A high DoI is required to maximize the energy density as well as to minimize the volume and cost of the system. The highest DoI is recorded in shaftless rotor layouts. The shaft covers only auxiliary functions and is not directly involved in the energy storage process. Furthermore, it raises the complexity, lowers the overall reliability and makes the dynamic performance of the flywheel worse by introducing critical resonances into the speed range of the system [22]. Therefore, straightforward advantages arise from a single-piece rotor in a shaftless configuration carrying out both the task of the energy accumulator and the rotating part of bearings and electric machine. Among the shaftless layouts, the outer-rotor configuration is the most favorable since all the mass is placed at the largest possible radius maximizing the rotational inertia with respect to the rotor weight. Furthermore, the stator components fit inside the rotor achieving a highly compact assembly [23]. Many studies in the literature are focused on this configuration. Kailasan et al. [24] proposed a layout where the stator lies in a hollow cylinder rotor made of two steel spline rings wrapped up with a light weight carbon-epoxy composite. The function of the steel spline rings is to connect the flywheel with the electric machine and the magnetic bearings. At IMS (TU Darmstadt) [23], a research team developed and tested a prototype of an outer-rotor flywheel where the rotor, made of high-strength fiber-reinforced composite material, is suspended on radial AMBs and a passive thrust magnetic bearing. These projects make use of a composite rotor because the match between the composite material and the outer-rotor shaftless layout maximizes the energy density of the system. As mentioned, the use of composite materials raises the cost of the system, but leads only to a small increment of the energy density since the rotor itself represents only a part of the overall weight [25]. Conversely, steel flywheels feature higher reliability, lower cost and easier manufacturing than composite designs. The energy density is not a primary concern in ground-based applications, which rely on the cost as well as on the system reliability [26,27]. This perspective has kept the research in shaftless steel flywheels alive. Amber Kinetics [25,28] proposed a shaftless flywheel with a solid steel rotor. This layout achieves a high DoI thanks to the rotating parts of the electric machine and bearings which contribute significantly to the rotating inertia. Li et al. [26] developed a shaftless flywheel with a solid disk steel rotor. The design is mass saving due to the shaftless layout and the magnetic levitation, which provides five-axes active magnetic suspension with a single bearing.

This paper provides a study of a shaftless outer-rotor FESS which uses a semi-hard magnetic steel rotor, with a five-axes active magnetic suspension and driven through two radial hysteresis self-bearing drives. In the context of the research on shaftless flywheels, the main contribution of this work is the study of the novel layout featuring a rotor made of semi-hard magnetic material (SHMM). Such material presents large mechanical strength allowing large peripheral speed and, consequently, high energetic performance of the FESS. Its magnetic features can be used to achieve suspension forces and, furthermore, the significant magnetic hysteresis that it exhibits can be exploited to yield mechanical torque [29,30]. This enables the use of the self-bearing hysteresis drives which offer the proposed layout low mechanical complexity, reduced material and manufacturing costs along with a lower component count. The flywheel is expected to store around 600 Wh at the maximum speed of 18,000 rpm. The mass of the rotor is 61.2 kg.

The authors intend to present the novel FESS idea focusing on the mechanical design, the assessment of the flywheel energetic performance, and the analysis of the radial magnetic levitation performance. The radial self-bearing drive is not investigated in detail; only its magnetic bearing operation is analyzed to assess the performance of the magnetic suspension based on the SHMM. For this purpose, the model of the suspended flywheel is generated. The electromagnets providing the levitating forces are modeled through a one-dimensional approach, where the Jiles–Atherton model is used to represent the magnetic hysteresis of the rotor material.

The paper is structured as follows. The proposed FESS layout is introduced in Section 2. Section 3 presents the flywheel mechanical design and analysis, composed of the energy and rotordynamic studies. The modeling of the radial magnetic suspension is presented in Section 4 and the results are analyzed in Section 5.

2. System Layout

The proposed layout of FESS is shown in Figure 1. It features a shaftless outer-rotor flywheel with a semi-hard magnetic steel rotor ①.

The flywheel is driven through two homopolar hysteresis self-bearing drives ② providing also the active magnetic suspension in the radial direction. The proposed architecture does not include any permanent magnet, unlike the common homopolar solution found in the literature [14,31]. The thrust AMB ③ is interfaced with two flanges ④ made of soft-magnetic material acting as an interface for the thrust backup bearings ⑤. The rotor is a hollow cylinder mounting two flanges. Such a smooth design is advantageous for flywheel manufacturing as well as its dynamic behavior, consequently requiring a smaller effort for the radial suspension. In this highly integrated layout, the flywheel acts as energy accumulator as well as rotating part of the self-bearing drives. The search for a more compact system is addressed by an outer-rotor solution, which allows all the stator components to fit into the hollow rotor. The drawback of the hollow shape for a steel flywheel is the reduced storage capability along with a lower energy density with respect to a solid steel rotor solution [32]. This is not a primary concern since the integrated layout and the self-bearing solution raise the energy density. Furthermore, the compactness of the system improves its specific energy thus making the layout advantageous especially for ground-based applications.

The rotor is made of FeCrCo 48/5, a semi-hard magnetic high-strength alloy. It is advantageous with respect to other SHMMs because of the lower content of cobalt (around 5%) along with a good combination of mechanical strength and energy density [29,33]. SHMMs feature a large magnetic hysteresis behavior that can be exploited to yield mechanical torque. This principle is used in hysteresis motors, whose application is, however, narrowed nowadays due to the limited power density with respect to other electric machine solutions. Nevertheless, hysteresis machines, thanks to the high-strength of SHMMs, are suitable for high and ultra-high speed applications thus allowing unmatchable power density. Furthermore, they do not require any additional components on the rotor since the rotating part of such machines is a simple ring made of semi-hard magnetic alloy. The hysteresis machine can also operate in the generator mode. The semi-hard magnetic material features a residual magnetization, which enables the machine to work as a generator. In this case, its operation is similar to a permanent-magnet generator.

The magnetic features of the FeCrCo 48/5 SHMM are also used to achieve levitating forces for the radial centering. Usually, the rotating part of magnetic bearings is made of soft-magnetic material. Although the use of SHMM could lead to high hysteresis losses as well as low force capability due to its large hysteretic behavior and poor magnetic permeability, these issues are not a primary concern for the proposed flywheel concept. The flywheel has to be vertically oriented, hence its weight is not supported by the radial suspension. Furthermore, the force due to residual unbalance is not significant because the system is designed to make the rotor work in a self-centered condition in the whole speed range. Therefore, low force is required to the radial suspension in normal operation. Moreover, the homopolar layout for the radial electric machine enables the reduction of iron losses.

The proposed layout for the radial self-bearing drive is shown in Figure 2. The motor-generator function exploits a three-phase distributed winding scheme on a twelve-slot stator (Figure 2a). The bearing function uses the homopolar configuration illustrated in Figure 2b. A single tooth for the motor- generator function in the transversal $xy$-plane is split in two teeth in each of the longitudinal $xz$-actuation plane for the bearing function. The magnetic fluxes responsible for the motor-generator and bearing functions lie on the $xy$-plane and on the longitudinal $xz$-actuation plane, respectively. Figure 2c shows the isometric view of the self-bearing drive, where the three-phase winding of the motor-generator function is entirely represented. Conversely, only the coils responsible for the force generation along one semi-axis are reported (each bearing actuation semi-axis involves three teeth of the $xy$-plane). A refined analysis is required to assess the performance of the proposed configuration. However, this goes beyond the aims of this paper.

3. Mechanical Design, Energy and Rotordynamic Analysis

The study of the flywheel is carried out analyzing the stress, energy and rotordynamic behavior of the system. The physical and mechanical properties of the FeCrCo 48/5 SHMM are reported in

Table 1. Properties of the FeCrCo 48/5 SHMM.

Symbol	Quantity	Value
$\nu$	Poisson ratio	0.3
$S_{UT}$	Tensile strength	1150 MPa
$\rho$	Mass density	7600 kg/m${}^3$

.

3.1. Mechanical Design and Energy Analysis

During its life, a flywheel undergoes a variety of load conditions. Among these, only the centrifugal loading is considered in this study since other load conditions lead to negligible effects [32]. The outer radius $r_o$ (150 mm) and the storage capability $E_{max}$ (600 Wh) are considered as design specifics. Therefore, the mechanical design has to be completed by defining the height h and the inner radius $r_i$.

Under the assumption of plane stress, the stress field of an axisymmetric body is a function of the radius r only. This assumption is valid for thin disks and offers a simple approach yielding a closed form solution. The plane stress field features only radial (${\sigma }_r$) and hoop (${\sigma }_h$) stress, whereas shear stress vanishes because of the problem axisymmetry. In the case of a hollow disk, the stress field is described by [32].

$$\begin{array}{c}{\sigma }_r(r,\beta )=\rho {\mathsf{\Omega }}^2{F}_r(r,\beta )\hfill \\ {\sigma }_h(r,\beta )=\rho {\mathsf{\Omega }}^2{F}_h(r,\beta ),\hfill \end{array}$$

(1)

where

$$\begin{array}{c}{F}_r(r,\beta )=\frac{3+\nu }{8}{r}_o^2\left(1+{\beta }^2\left(1-{\left(\frac{r_o}{r}\right)}^2\right)-{\left(\frac{r}{r_o}\right)}^2\right)\hfill \\ F_h(r,\beta )=\frac{3+\nu }{8}{r}_o^2\left(1+{\beta }^2\left(1+{\left(\frac{r_o}{r}\right)}^2\right)-\frac{1+3\nu }{3+\nu }{\left(\frac{r}{r_o}\right)}^2\right),\hfill \end{array}$$

(2)

$\mathsf{\Omega }$ is the rotational speed, $\nu$ and $\rho$ are the Poisson ratio and the density of the rotor material, respectively, and $\beta$ is the ratio between the inner and the outer radius $r_i/r_o$. Equations (1) and (2) hold only for $\beta \in (0;1)$. By considering the Von Mises criterion reduced to the principal plane stress case, the equivalent stress ${\sigma }_{eq}$ is

$${\sigma }_{eq}(r,\beta )=\rho {\mathsf{\Omega }}^2\overline{F}(r,\beta ),$$

(3)

where

$$\overline{F}(r,\beta )=\sqrt{F_r{(r,\beta )}^2-F_r(r,\beta )F_h(r,\beta )+F_h{(r,\beta )}^2}.$$

(4)

The maximum equivalent stress occurs at the inner radius $r_i$. The burst speed happens when the maximum equivalent stress matches the tensile strength of the rotor material. From Equation (3), by equaling the equivalent stress at the inner radius to the tensile strength $S_{UT}$

$${\mathsf{\Omega }}_{Burst}\left(\beta \right)=\sqrt{\frac{S_{UT}}{\rho \overline{F}(r_i,\beta )}}.$$

(5)

Since ${\mathsf{\Omega }}_{Burst}$ is related to the burst of the flywheel, the maximum speed in normal operation needs to consider a safety factor $S{F}_{\mathsf{\Omega }}={\mathsf{\Omega }}_{Burst}/{\mathsf{\Omega }}_{Max}$. One acceptability criterion used for turbomachinery considers that the rotating components must be able to support a speed up to $40\%$ greater than the maximum speed in normal operation. Hence a safety factor $S{F}_{\mathsf{\Omega }}$ of $1.4$ is considered as typical [32].

The energy storage capability of the flywheel, which is required to be $600$ Wh, is equal to its kinetic energy at the maximum speed

$$E_{Max}=\frac{1}{2}{I}_p\left(\beta \right){\left(\frac{{\mathsf{\Omega }}_{Burst}\left(\beta \right)}{S{F}_{\mathsf{\Omega }}}\right)}^2,$$

(6)

where

$$I_p\left(\beta \right)=\frac{1}{2}\pi h\rho r_o^4(1-{\beta }^4)$$

(7)

is the polar moment of inertia of the flywheel. The height h that gives a storage capability $E_{max}$ at the different values of the ratio $\beta$ is obtained by substituting Equations (5) and (7) into Equation (6)

$$h\left(\beta \right)=\frac{4E_{Max}}{\pi r_o^4(1-{\beta }^4)}\frac{\overline{F}(r_i,\beta )}{{\sigma }_{lim}},$$

(8)

where the ratio ${\sigma }_{lim}=S_{UT}/S{F}_{\mathsf{\Omega }}^2$ is the maximum equivalent stress occurring in normal operation. For determining the ratio $\beta$ and, consequently, the height h through Equation (8), a constraint on the ratio between the polar and the transversal moment of inertia has been considered, where the transversal moment of inertia ($I_t$) with respect to the barycentric axis is

$$I_t\left(\beta \right)=\frac{1}{4}\pi \rho r_o^{2}h\left(\beta \right)\left[\frac{h{\left(\beta \right)}^2}{12}+r_o^2(1+{\beta }^2)\right](1-{\beta }^2)$$

(9)

The ratio $I_p/I_t$ covers a relevant role in the unbalance response of the spinning flywheel since it determines the occurrence of the conical critical speed. Values sufficiently lower than 1 are usually preferred to reduce the conical critical speed, thus widening the speed range in which the rotor is self-centered (condition that reduces the request to the radial bearings) [34,35]. Figure 3 shows both the height h and the ratio $I_p/I_t$ versus the ratio $\beta$. The selected ratio $I_p/I_t$ is 0.5, corresponding to $\beta$ equal to 0.9. The resulting height h is 600 mm. The designed flywheel dimensions are summarized in

Table 2. Dimensions of the flywheel.

Symbol	Quantity	Value
$r_i$	Inner radius	135 mm
$r_o$	Outer radius	150 mm
h	Height	600 mm

.

A maximum spin speed of 18,000 rpm is obtained by substituting the achieved value of $\beta$ into Equation (5), divided by the safety factor $S{F}_{\mathsf{\Omega }}$. The minimum speed is selected to avoid excessive current values as well as high frequency drop [36]. It is obtained by setting the Depth of Discharge (DoD) to 89%, which is commonly used in the literature [37]. Therefore, the working speed range of the flywheel is 6000–18,000 rpm.

The energy density, which is the ratio between the storage capability and the mass m of the flywheel, can be written

$$\frac{E_{Max}}{m}=K\frac{{\sigma }_{lim}}{\rho },$$

(10)

where

$$m=\pi h\rho r_o^2(1-{\beta }^2),K=\frac{1}{4}{r}_o^2\frac{1+{\beta }^2}{\overline{F}(r_i,\beta )}.$$

(11)

The shape factor K quantifies the shape efficiency of the flywheel. It depends only on the rotor geometry and on the employed failure criterion. Thus, it can be calculated that the proposed flywheel features a mass of 61.2 kg and a mass-based energy density of 9.8 Wh/kg. The achieved shape factor is 0.468. This value could be increased toward the maximum value for hollow disk flywheels of 0.5 by raising the ratio $\beta$ toward the ideal value of 1. In the case of solid disk flywheel, the shape factor is equal to 0.606 [32]. This explains why the solid cylinder shape for a steel flywheel leads to higher energetic performance with respect to the hollow cylinder shape, as remarked in Section 2.

As mentioned above, the plane stress assumption leads to valid results for the case of thin disk. Therefore, the long cylinder case requires a more accurate analysis to evaluate the three-dimensional stress field and assess the value of the plane stress approach. The analysis is carried out by means of a finite element (FE) model of the designed flywheel, which is built in ANSYS^® Mechanical APDL 15.0 (Pennsylvania, Canonsburg, USA). Due to the axisymmetry of both geometry and loads, a two-dimensional model is sufficient. The axisymmetric formulation of the PLANE82 element is used to reduce the size of the problem. The rotor, spun at 18,000 rpm, is constrained only in the axial direction, since other translations are constrained by the axial symmetric formulation of the problem. Resulting stresses from the plane stress and the FE approaches are evaluated along the radius and compared in Figure 4. The difference between the two solutions is lower than the $0.5\%$. This comparison demonstrates that the analytical solution well describes the stress field under the centrifugal loading.

3.2. Rotordynamics

Since the rotor has a low radial thickness, vibrations of the cross-section of the rotating cylinder may occur at low frequency. Extensive discussions on the dynamics of thin-walled rotating cylinders can be found in the literature [38,39,40].

A three-dimensional modal analysis of the rotating flywheel was carried out by means of a model built in ANSYS^® Mechanical APDL 15.0. Even if the rotor is axisymmetric, the three-dimensional model is needed to evaluate the three-dimensional mode shapes. ANSYS Solid95 elements (twenty-node brick elements) are used because they embed both gyroscopic and Coriolis effects. The latter fulfills an important role in the dynamics of thin-walled rotating cylinders.

The goal of the study is to check if any of the deformation modes of the flywheel affect its speed range. The free–free boundary condition is chosen. It is a conservative condition, since any constraint would stiffen the system thus increasing the magnitude of the natural frequencies. The natural frequencies $\omega$ at different spin speed $\mathsf{\Omega }$ are computed and reported in the Campbell diagram (Figure 5). The pure radial and the radial shearing mode shapes are shown in Figure 6. They involve vibrations of the cross-section of the rotor. The Campbell diagram shows that the radial shearing mode is significantly affected by the Coriolis effect. It leads to the bifurcation of the natural frequencies by lowering the forward (FW) frequency and raising the backward (BW) one. It is important to notice that the Campbell diagram of the transversal cross-section modes has a different behavior compared to the one related to the classical flexural mode shapes [38]. The centrifugal effect tends to raise both the FW and the BW frequencies of the cross-section modes. Therefore, it has been neglected in this analysis. As such, the forward natural frequency of the radial shearing mode decreases with increasing the spin speed $\mathsf{\Omega }$. Nevertheless, it does not cross the synchronous $\omega =\mathsf{\Omega }$ line and, consequently, neither resonances nor instability phenomena occur. The pure radial mode is affected slightly by the Coriolis effect. FW and BW frequency bifurcation is not pronounced, their magnitude difference is negligible and both frequencies are almost constant in the considered speed range. The frequency of the first deformation mode (the pure radial) is nearly twice the top speed of the flywheel. Therefore, the rotor can be considered as a rigid body in the study of the magnetic suspension.

4. Radial Active Magnetic Suspension

Since the goal is the performance assessment of the radial magnetic suspension based on the SHMM, only the magnetic bearing operation of the self-bearing drive is analyzed. For this purpose, the study considers two homopolar radial AMBs instead of the radial self-bearing drives. This enables to follow the common modeling approach used for AMBs [3]. Since an open-loop active magnetic suspension is inherently unstable, a Proportional-Derivative (PD) control law is adopted to stabilize the system.

Table 3. Main parameters of the radial AMB.

Symbol	Quantity	Value
$A_p$	Pole area	2010 mm${}^2$
$g_0$	Nominal magnetic gap	0.8 mm
$l_r$	Length of the magnetic flux path in the rotor	36 mm
N	Turns per electromagnet	174
$i_0$	Bias current	2.5 A
$i_{max}$	Maximum continuous current	5 A

provides the main parameters of the considered radial AMBs.

Modeling

The present section describes the modeling of the magnetically suspended flywheel. Figure 7 shows a schematic of the rigid rotor on AMBs.

The dynamics of the rigid rotor on AMBs at constant spin speed [3,34] is governed by

$$\mathbf{M}\ddot{\mathbf{q}}+\mathsf{\Omega }\mathbf{G}\dot{\mathbf{q}}={\mathbf{f}}_{\mathbf{unb}}+{\mathbf{f}}_{\mathbf{ext}},$$

(12)

where

$$\mathbf{M}=diag\left\{\begin{array}{c}{I}_t\\ m\\ I_t\\ m\end{array}\right\},\mathbf{G}=\left[\begin{array}{cccc}0& 0& I_p& 0\\ 0& 0& 0& 0\\ -I_p& 0& 0& 0\\ 0& 0& 0& 0\end{array}\right],{\mathbf{f}}_{\mathbf{unb}}={\mathsf{\Omega }}^2\left\{\begin{array}{c}-\chi (I_p-I_t)cos\mathsf{\Omega }t\\ mϵcos\mathsf{\Omega }t\\ -\chi (I_p-I_t)sin\mathsf{\Omega }t\\ mϵsin\mathsf{\Omega }t\end{array}\right\},$$

m, $I_t$ and $I_p$ are the rotor mass, transversal moment of inertia and polar moment of inertia, respectively, and t is the simulation time. $mϵ$ and $\chi (I_p-I_t)$ are the static and couple unbalance, which are computed in accordance with the standards reported in [34] (the quality grade suggested for flywheels is G 6.3). The values of these parameters are reported in

Table 4. Mechanical parameters of the rotor.

Symbol	Quantity	Value
m	Mass	61.2 kg
$I_t$	Transversal moment of inertia	2.46 kgm${}^2$
$I_p$	Polar moment of inertia	1.25 kgm${}^2$
$\chi (I_p-I_t)$	Couple unbalance	4.50 × 10−5 kgm2
mϵ	Static unbalance	2.05 × 10−4 kgm

.

The vector $\mathbf{q}={\left[\eta x\psi y\right]}^T$ contains the generalized coordinates of the rotordynamics (x and y are the translational degrees of freedom, $\eta$ and $\psi$ are the rotational ones, as illustrated in Figure 7), and ${\mathbf{f}}_{\mathbf{ext}}$ is the vector of the external forces acting on the rotor:

$$\begin{array}{c}{\mathbf{f}}_{\mathbf{ext}}={\mathbf{f}}_{\mathbf{dist}}+{\mathbf{B}}_{\mathbf{T}}{\mathbf{f}}_{\mathbf{AMB}},\\ {\mathbf{B}}_{\mathbf{T}}=\left[\begin{array}{cccc}\mathrm{a}& \mathrm{b}& 0& 0\\ 1& 1& 0& 0\\ 0& 0& \mathrm{a}& \mathrm{b}\\ 0& 0& 1& 1\end{array}\right],{\mathbf{f}}_{\mathbf{AMB}}=\left\{\begin{array}{c}{f}_{Ax}\\ f_{Bx}\\ f_{Ay}\\ f_{By}\end{array}\right\},\end{array}$$

(13)

where $\mathrm{a}=-220\mathrm{mm}$ and $\mathrm{b}=220\mathrm{mm}$ are the coordinates of the radial AMBs with respect to the origin of the $xyz$-reference frame. ${\mathbf{f}}_{\mathbf{dist}}$ is the vector of the external disturbances acting on the flywheel. The vector ${\mathbf{f}}_{\mathbf{AMB}}$ contains the forces provided by the AMBs, where $f_{ij}$ is the force of the AMB i along the j-axis. Two counteracting electromagnets, driven in the differential driving mode [3], are employed along each axis of the AMB. By considering the j-axis of the AMB i, $f_{e,ij}^{+}$ and $f_{e,ij}^{-}$ are defined as the forces provided by the electromagnets supplied with the currents $i_0+i_{c,ij}$ and $i_0-i_{c,ij}$, respectively, where $i_0$ is the bias current and $i_{c,ij}$ is the control current at the j-axis of the AMB i. The latter is determined by the PD controller, which receives as input the displacement of the flywheel at the AMB location [3]. The parameters of the PD controller are reported in

Table 5. Parameters of the PD controller in the form of $C\left(s\right)=P+\frac{D}{1+{\tau }_{d}s}$.

Symbol	Quantity	Value
P	Proportional gain	12,039 $\mathrm{A}/\mathrm{m}$
D	Derivative gain	$111\mathrm{A}\mathrm{s}/\mathrm{m}$
${\tau }_d$	Time constant	$3.3\mathrm{m}\mathrm{s}$

.

The force $f_{ij}$ is computed as

$$f_{ij}=f_{e,ij}^{+}-f_{e,ij}^{-}$$

(14)

A one-dimensional approach is used to achieve the force $f_e$ produced by the electromagnet. The magnetic circuit of the electromagnet is sketched in Figure 8, where $l_r$ is the length of the magnetic flux path in the rotor and g is the magnetic gap. The current i flows into the coil and forces the magnetic flux $\mathsf{\Phi }$, which is constant if flux leakages are neglected, in the magnetic circuit.

Ampère’s law can be applied to the circuit:

$$H_s{l}_s+2H_{g}g+H_r{l}_r=Ni,$$

(15)

where $l_s$ is the length of the flux path in the stator, N is the number of turns of the coil and H is the magnetic field strength (the subscript s refers to the stator, g to the magnetic gap, r to the rotor). The core of the electromagnet is made of laminated silicon steel, whose magnetic behavior can be approximated linearly. By introducing the linear magnetic characteristics of soft-iron and free-space, Equation (15) gives

$$\frac{\mathsf{\Phi }}{A_s}\frac{l_s}{{\mu }_0{\mu }_s}+2\frac{\mathsf{\Phi }}{A_p}\frac{g}{{\mu }_0}+H_r{l}_r=Ni,$$

(16)

where $A_s$ and $A_p$ are the cross-sectional areas of stator and magnetic gap, respectively, ${\mu }_0$ is the permeability of free-space and ${\mu }_s$ is the relative permeability of silicon steel. The latter is significantly larger than the permeability of the free space, therefore the term referring to the stator can be neglected. If the cross-sectional areas along the flux path do not change, the magnetic flux density B is constant ($B_s$ = $B_a$ = $B_r$ = B) and Equation (16) gives

$$B=\frac{{\mu }_0}{2g}Ni-\frac{{\mu }_0{l}_r}{2g}{H}_r.$$

(17)

It is the load curve of the rotor material, which is a negative slope line in the B-H diagram whose distance from the origin mainly depends on the current supplied to the coil. The operating point of the rotor material yields from the intersection between the load curve and the magnetic characteristic of the rotor material. The procedure is similar to that adopted for the operating point of a permanent magnet [41]. In the present case, the demagnetization line of the permanent magnet is substituted by the magnetic characteristic of the rotor material. The latter exhibits a magnetic hysteretic behavior, which is represented through the Jiles–Atherton model [42,43]. The magnetic constitutive relation of the rotor material is

$$\begin{array}{c}\hfill B={\mu }_0(H_r+M_r),\end{array}$$

(18)

where the magnetization $M_r$ is obtained from the first order differential equation of the Jiles–Atherton model

$$\begin{array}{c}\frac{\mathrm{d}{M}_r}{\mathrm{d}{H}_r}=\left[c\frac{\mathrm{d}{M}_{r,an}}{\mathrm{d}{H}_{r,eff}}+\left(1-c\right)\frac{{\delta }^{\star }}{k\delta }\left(M_{r,an}-M_{r,irr}\right)\right]\frac{\mathrm{d}{H}_{r,eff}}{\mathrm{d}{H}_r},\\ \delta =sgn\left(\mathrm{d}{H}_{r,eff}\right),{\delta }^{\star }=\left\{\begin{array}{cc}1\hfill & if\left(M_{r,an}-M_{r,irr}\right)\mathrm{d}{H}_{r,eff}>0\hfill \\ 0\hfill & \mathrm{otherwise}\hfill \end{array}\right..\end{array}$$

(19)

The parameter c describes the magnetization reversibility and k is the pinning coefficient. $M_{r,an}$ is the anhysteretic magnetization, and $H_{r,eff}$ is the effective magnetic field

$$\begin{array}{c}\hfill H_{r,eff}=H_r+\alpha M_r,\end{array}$$

(20)

where $\alpha$ represents the interdomain coupling. The anhysteretic magnetization $M_{r,an}$ is expressed through the Langevin function

$$\begin{array}{c}\hfill M_{r,an}=M_s\left[coth\left(\frac{H_{r,eff}}{a}\right)-\frac{a}{H_{r,eff}}\right].\end{array}$$

(21)

The parameters needed to represent magnetic characteristic of the FeCrCo 48/5 (

Table 6. Parameters of the Jiles–Atherton model for the FeCrCo 48/5 SHMM.

Symbol	Quantity	Value
$\alpha$	Local field factor	$0.15$
c	Domain rotation loss	$0.20$
a	Langevin parameter	$9.94\times {10}^4\mathrm{A}/\mathrm{m}$
k	Pinning	$5.16\times {10}^4\mathrm{A}/\mathrm{m}$
$M_S$	Saturation magnetization	$2.02\times {10}^6\mathrm{A}/\mathrm{m}$

) have been identified by means of the least-squares optimization method proposed by Kis and Iványi [44]. The identification has been based on the material hysteresis loop measured on FeCrCo 48/5 SHMM samples by means of a vibrating-sample magnetometer (VSM). Figure 9 shows the comparison between the experimental and the numerical hysteresis loop. Equation (19) requires the initial condition of the magnetization $M_r$, which is computed through the constitutive relation (Equation (18)) by providing the initial magnetic field $H_{r,0}$ and the initial magnetic field density $B_0$. Therefore, the magnetic characteristic of the rotor material can be written as

$$B=B(H_r,H_{r,0},B_0).$$

(22)

The intersection between Equations (17) and (22) is solved numerically at each step time of the simulation and gives the operating point of the rotor material (${\overline{H}}_r,\overline{B}$). Note that the numerical implementation of the model requires that the initial condition of Equation (22) at the step time n is the operating point of the rotor material at the preceding step time (${\overline{H}}_{r,n-1},{\overline{B}}_{n-1}$). The force that the electromagnet generates yields from the Maxwell stress tensor [41]. The one-dimensional approach considers a constant flux density across the magnetic gap, hence the Maxwell stress tensor gives

$$f_e=\frac{A_p}{{\mu }_0}{\overline{B}}^2.$$

(23)

The force $f_e$ is computed for both the counteracting electromagnets and then it is introduced into Equation (14) to achieve the force that the AMB generates along the considered direction.

In the case of the spinning flywheel, the model needs to be modified to consider the residual magnetization effect of the hysteresis phenomenon. Figure 10 displays a sketch of a generic rotor with the four electromagnets of the AMB. A small cube of material is highlighted. At the time $t_n$, it interacts with Electromagnet 1 and, consequently, the material operating point is ${({\overline{H}}_{r,n},{\overline{B}}_n)}_1$. When the rotor rotates, the cube of material leaves the electromagnet 1 but keeps its magnetization. As soon as it approaches the region of Electromagnet 2, its operating point will depend on the operating point at the previous electromagnet. The modeling approach considers that the step time $dt$ is one quarter of the rotational period. Therefore, it is the time that the cube of material needs to move from one electromagnet to the next. Furthermore, the operating point of the rotor material is assumed not to change in between the two electromagnets. This is implemented by providing the proper initial condition at Equation (22). By considering Electromagnet k at the time $t_{n+1}$, the initial condition to be provided is the operating point ${({\overline{H}}_{r,n},{\overline{B}}_n)}_{k-1}$ of the preceding Electromagnet $k-1$ at the previous step time $t_n$.

5. Results

The model of the suspended flywheel is used to simulate the response to an external disturbance at $\mathsf{\Omega }$ = 0 rpm and the response to the unbalance at the maximum speed $\mathsf{\Omega }$ = 18,000 rpm. The former allows investigation of the magnetic hysteresis effect on the magnetic suspension. The latter is needed to assess the performance of the system at the maximum speed of the flywheel.

5.1. Response to the External Disturbance at Standstill

The system at standstill is disturbed by exerting a force step to the rotor center of mass. Therefore, the load is divided symmetrically on the two actuation stages. The response of only one AMB along the perturbed direction is reported since both the AMBs respond in the same way. The response of the system is compared with the response of a similar system, which approximates linearly the FeCrCo 48/5 magnetic behavior. For this purpose, the magnetic characteristic of the rotor material (Equation (22)) is substituted by a linear relationship $B={\mu }_0{\mu }_r{H}_r$, where ${\mu }_r=13$ is yielded from the slope of the first magnetization curve of the FeCrCo 48/5.

Figure 11a,b show the resulting displacement of the rotor at the location of the radial AMBs for step amplitudes of 100 N and 300 N, respectively. One of the major effects of the hysteresis phenomenon is that the AMB is less responsive to a current change, which can be attributed to a lag between the input and the output of the AMB. This effect is evident from the response of the system; the force generated by the AMB (Figure 11c,d rises more slowly with respect to the system with linearized rotor material whereas the current change is quicker (Figure 12). The current in the case of the 300-N step amplitude saturates at 10 A. Another important effect of the hysteresis is the residual magnetization of the rotor material. This aspect is investigated at various step amplitudes by considering the steady-state response of the force provided by the most loaded electromagnet (Figure 13a). Such force, in the case of hysteretic rotor material, is composed of two contributions: the residual force (dash-dotted line), which is due to the residual magnetization of the rotor material, and the current-related force (dashed line), which is due to the supplied current. The residual magnetization contributes mainly to the total force (solid line). This force is higher than that required by the load (half of the load applied to the center of mass) because it includes the balancing of the counteracting electromagnet force due to the residual magnetization. On the other hand, in the linearized rotor material case (solid line with circles) this effect is not present and, consequently, the value of force is lower. The current requested by the electromagnet (Figure 13b) in the hysteretic rotor material case (solid line) is significantly lower than the current needed in the linear rotor material case (solid line with circles). The residual magnetization and hence the residual force increases with the amplitude of the transient current overshoot (dash-dotted line), whose value is reported in Figure 13b. When the step amplitude is lower than 150 N, the increase in the total force is compensated by the residual force. The resulting control current ($i_c$) that has to be summed to the bias ($i_0$) is then small. As the current overshoot starts saturating, the force requested to the current is higher. It is worth mentioning that, as stated in [45], the SHMM-based magnetic suspension features a reduced force capability when compared to common solutions with soft iron.

5.2. Response to the Unbalance at Maximum Speed

The response to the unbalance is evaluated at 18,000 rpm. The material of the rotor is considered magnetically virgin when the simulation starts, hence the homogeneous condition $(H_0,B_0)=(0,0)$ is provided to each electromagnet. The initial conditions of the generalized coordinates of the lateral dynamics ($\mathbf{q}$) are homogeneous. The unbalance force appears suddenly as soon as the simulation starts, hence the rotor is subject to a step force and a transient response is present. The force step at the starting time needs to be avoided in order to limit the transient time and achieve the steady state response in a shorter simulation time. For this purpose, the unbalance signal is multiplied by a gain which ramps linearly from 0 to 1 in 0.02 s.

Figure 14a,b show the maximum displacement of the rotor at the location of the radial AMBs. The maximum amplitude at steady state occurs at the AMB B. It is one order of magnitude lower than the nominal magnetic gap, hence the response is satisfactory. The forces provided by each AMB are reported in Figure 14c,d. Figure 14e,f show the currents required for the generation of such forces. Their values are smaller than the allowable continuous current that can be supplied to the electromagnets (Table 3).

6. Conclusions

This paper presented a novel FESS layout featuring a shaftless flywheel made of FeCrCo 48/5 SHMM on self-bearing hysteresis drives. The study of the stress due to the centrifugal effect led to the design as well as an assessment of the flywheel maximum speed and its energetic features. The rotordynamic analysis, performed by means of a three-dimensional FE approach, showed that the rotor deformation modes were well above the maximum speed of the flywheel, which was thus considered as a rigid body in the consecutive analyses. The performance of the radial SHMM-based magnetic suspension was investigated by considering two homopolar radial AMBs instead of the radial self-bearing drives. The modeling of the suspended flywheel was presented and the responses to both the external step force at standstill and to the unbalance at the maximum speed were achieved. These analyses allowed to highlight the main effects of the hysteresis of the rotor material on the magnetic suspension and to prove its satisfactory performance at the maximum speed of the flywheel.

The presented study led to the evaluation of the flywheel energetic features and the assessment of the radial magnetic suspension performance. Future work can include a detailed self-bearing drive study in order to consider both its bearing and electric machine operations.