Analytical Approximation of the Stress Function for Conical Flywheels

1. Introduction

Calculating the stress of variable-thickness flywheels using an equation is of absolute importance in mechanical engineering. Most flywheel cross sections are designed and optimized using FEA software as a consequence of the lack of an equation to determine the stresses.

Flywheels have multiple applications, such as space technology, renewable energy and transportation systems, and power networks, among others [1,2,3,4]. Flywheels are the main component of flywheel energy storage systems (FESS), an electromechanical device that converts kinetic rotational energy into electrical energy [2]. FESS are widely used for their flexibility to store energy and later convert and manipulate it [2]. They offer high efficiency, large amounts of instantaneous power, fast response, low maintenance costs, long-lasting services, in addition to multiple environmental benefits [1,3,4].

1.1. Background

In analyzing and designing an FESS, optimizing the size, shape, and topology of the flywheel is of great value, since these parameters are directly related to the energy storage capacity of the FESS [2]. There are many researchers who investigate the analysis and optimization of flywheel size, shape, and topology [5,6,7]. On the one hand, most of them use numerical methods, such as the finite element method (FEM), or finite element analysis (FEA). Kailisan [4], for example, performed an FEA on a composite flywheel that integrates carbon fiber and a steel spline ring to limit the tangential stress that would allow high-speed operation. The authors show positive results, including the fact that, because of the steel rings, centrifugal growth would not induce the rotor parts to separate. A shape optimization FEM was implemented by Jiang [6]. The method relies on the parametric geometry modeling method and the downhill simplex method. The results show that the energy density and working safety performance could be significantly improved by implementing this shape optimization method. Furthermore, Jiang [7] implemented the variable density method and FEM to optimize the topology of a flywheel for energy density. The results show a 56.7% increase in the energy density of a constant thickness flywheel. Gao [8] presents a novel FESS, with the flywheel design being the main novelty. It consists of a double hub and rim that yields a mass reduction and an increase in inertia. The results show no separation between the rim and the dual hubs at high speeds. Kale [5] introduces an optimization method that eliminates the constraint of arbitrarily large volumes, used in traditional kinetic energy methods. The method, based on the maximization of specific energy, achieved a 15.8% increase in specific energy, compared to the kinetic energy method.

On the other hand, there are researches that have proposed analytical solutions for variable thickness flywheels. But this is not an easy task. Wen [9], for example, analyzed the stresses in an anisotropic flywheel based on plane stress. After the radial and tangential stresses were obtained, the location of the maximum radial stress was derived using the extreme point method and calculated using the Newton method. The results show a fine approximation which could be used in the design of composite flywheels. To the best of our knowledge, at the time of developing this work, the equilibrium equation for a flywheel with a linearly variable thickness had not been solved. Although Semsri [2] uses the solution provided by Ugural [10] to analyze a conical flywheel, Ugural’s solution is clearly and explicitly provided for a hyperbolic profile.

In the next section, the equilibrium equation is solved, assuming that the flywheel’s thickness varies linearly. A compact and easy-to-handle equation is given and further analyzed. The equation obtained was used to compute the stresses developed in the flywheel, and these results were compared with the FEA results. In Section 3, the results of Section 2 are discussed.

2. Methodology

The design of a constant thickness flywheel is based mainly on two stresses: a radial stress, ${\sigma }_r$, and a tangential one, ${\sigma }_t$. ${\sigma }_r$ is given by [2]

$$\begin{array}{c}\hfill {\sigma }_r=\left(r_0^2+r_i^2-{\displaystyle \frac{r_0^2{r}_i^2}{r^2}}-r^2\right){\displaystyle \frac{3+\nu }{8}}\rho {\omega }^2\end{array}$$

(1)

and ${\sigma }_t$ by [2]

$$\begin{array}{c}\hfill {\sigma }_t=\left(r_0^2+r_i^2+{\displaystyle \frac{r_0^2{r}_i^2}{r^2}}-{\displaystyle \frac{1+3\nu }{3+\nu }}{r}^2\right){\displaystyle \frac{3+\nu }{8}}\rho {\omega }^2\end{array}$$

(2)

where

$r_o$ is the outer radius of the flywheel,

$r_i$ is the inner radius of the flywheel,

r is the radial distance at which the stress is to be computed,

$\nu$ is the flywheel’s material’s Poisson ratio,

$\rho$ is the density of the flywheel material and

$\omega$ is the flywheel’s angular velocity.

In Figure 1 the results of an FEM analysis, performed in a commercial FEA software, are shown. It can be seen that the results of the simulation converge with those computed analytically with Eq. 1 and 2. Figure 1a shows ${\sigma }_r$, Figure 1b shows ${\sigma }_t$, and Figure 1c shows the comparison between the simulation results and Eq. 1 and 2. Figure 1c shows the behavior of ${\sigma }_r$ and ${\sigma }_t$, using data from Table 1. It can be seen that the total stress, ${\sigma }_T=\sqrt{{\sigma }_r^2+{\sigma }_t^2}$, is much more influenced by ${\sigma }_t$. The maximum stress is about 14 kPa.

Although the stresses for a constant-thickness flywheel are relatively easy to find, the stresses in a variable-thickness flywheel are not. In this study, the main goal was to find an expression that would help to calculate the stresses in a variable-thickness flywheel. This is not an easy task, ever since there is no equation to compute the stresses in most variable thickness shapes. To find such an equation, the analysis started from the stress function [11] given by

$$\begin{array}{c}\hfill r^2{\displaystyle \frac{d^{2}F}{d{r}^2}}+r{\displaystyle \frac{dF}{dr}}-F+\left(3+\nu \right)\rho {\omega }^{2}t{r}^3-{\displaystyle \frac{r}{t}}{\displaystyle \frac{dt}{dr}}\left(r{\displaystyle \frac{dF}{dr}}-\nu F\right)=0\end{array}$$

(3)

Equation 3 is clearly not linear when t is described by a polynomial, as in the case of a conical flywheel. Therefore, it is necesary to solve Eq. 3 by other means. In Eq. 3, the fourth term derives from the rotational body force. That is, when such a term is considered, Eq. 3 is not homogeneous, while if it is not, then Eq. 3 is homogeneous.

The most simple case for a variable-thickness flywheel is when the thickness varies linearly along the radius. That being said, lets consider t given by

$$\begin{array}{c}\hfill t=A+Br\end{array}$$

(4)

Now, if the variation of t is small, that is, if $dt/dr$ is small, then the last term in Eq. 3 can be neglected and it becomes:

$$\begin{array}{c}\hfill r^2{\displaystyle \frac{d^{2}F}{d{r}^2}}+r{\displaystyle \frac{dF}{dr}}-F+\left(3+\nu \right)\rho {\omega }^{2}t{r}^3=0\end{array}$$

(5)

Since the last term in Eq. 5 yields to a fourth degree polynomial, considering Eq. 4, the particular solution to Eq. 5, $F_p$, should have the form

$$\begin{array}{c}\hfill F_p=C+Dr+E{r}^2+G{r}^3+H{r}^4\end{array}$$

(6)

Upon substitution of Eq. 6 and 4 in Eq. 5, $C=D=E=0$, $G=-A\left(3+\nu \right)\rho {\omega }^2/8$, and $H=-B\left(3+\nu \right)\rho {\omega }^2/15$. Therefore, $F_p$ is

$$\begin{array}{c}\hfill F_p=-\left(3+\nu \right)\left({\displaystyle \frac{A}{8}}+{\displaystyle \frac{Br}{15}}\right)\rho {\omega }^2{r}^3\end{array}$$

As for the homogeneous solution, $F_h$, it is noted that the homogeneous form of Eq. 5 is the same as the homogeneous one for a constant thickness flywheel [11]. Therefore, $F_h$ is

$$\begin{array}{c}\hfill F_h=c_{0}r+c_1/r\end{array}$$

(7)

where $c_0$ and $c_1$ are determined from the boundary conditions. In this case, ${\sigma }_r=0$ at $r=r_i$ and $r=r_o$.

Once $F_p$ and $F_h$ are determined, the general solution is given by

$$\begin{array}{c}\hfill F=F_h+F_p=c_{0}r+c_1/r-\left(3+\nu \right)\left({\displaystyle \frac{A}{8}}+{\displaystyle \frac{Br}{15}}\right)\rho {\omega }^2{r}^3\end{array}$$

(8)

and ${\sigma }_r$ and ${\sigma }_t$ can be obtained as [11]

$$\begin{array}{c}\hfill tr{\sigma }_r=F\end{array}$$

and

$$\begin{array}{c}\hfill t{\sigma }_t={\displaystyle \frac{dF}{dr}}+t\rho {\omega }^2{r}^2\end{array}$$

Upon applying boundary conditions, $c_0$ and $c_1$ are given by

$$\begin{array}{c}\hfill c_0={\displaystyle \frac{\left[\left(r_i+r_o\right)\left(r_i^3+r_o^3\right)+r_i^2{r}_o^2\right]8B+\left(r_i+r_o\right)\left(r_i^2+r_o^2\right)15A}{120\left(r_i+r_o\right)}}Z\end{array}$$

(9)

and

$$\begin{array}{c}\hfill c_1=-{\displaystyle \frac{8B{r}_i{r}_o\left(1+\frac{r_i}{r_o}+\frac{r_o}{r_i}\right)+15A\left(r_i+r_o\right)}{120\left(r_i+r_o\right)}}{r}_i^2{r}_o^{2}Z\end{array}$$

(10)

with $Z=\left(3+\nu \right)\rho {\omega }^2$.

If the flywheel’s thickness is considered to vary linearly along the radius, as defined by Eq. 4, then A and B can be defined as functions of $r_i$, $r_o$, internal thickness, $t_i$, and outer thickness, $t_o$, as

$$\begin{array}{c}\hfill B={\displaystyle \frac{t_o-ti}{2\left(r_o-ri\right)}}\end{array}$$

(11)

and

$$\begin{array}{c}\hfill A=t_i-B{r}_i\end{array}$$

(12)

With Eq. 8 and Eq. 9 through 12, the relation between the stresses in the flywheel and its geometry is explicit.

To validate Eq. 8, an FEM analysis was performed in an FEM software and the results were compared with those obtained with Eq. 8 and are presented in Figure 2 and Figure 3. The results of Figure 2 were obtained using the parameters given in Table 2. From Table 2 it is evident that the flywheel’s thickness variation is too small. Therefore, an analysis was performed to determine the interval of values for the thickness relation $t_{o,i}=t_o/t_i$, for which Eq. 8 is valid. In other words, for what values of $t_{o,i}$ does Eq. 8 give an acceptable approximation of ${\sigma }_r$ and ${\sigma }_t$ in a variable-thickness flywheel with t given by Eq. 4.

Both ${\sigma }_r$ and ${\sigma }_t$ were obtained with FEA along two lines, line 1 and line 2, depicted in Figure 2a, for $1\leqq t_{o,i}\leqq 4.5$. Furthermore, Figure 2b shows the Pearson correlation factor, P, between ${\sigma }_r$ and ${\sigma }_t$ obtained with FEA and Eq. 8. The figure shows a high correlation for both stresses in both lines, between Eq. 8 and FEA; Figure 2c shows the mean absolute percentage error (MAPE). It shows that the MAPE behaves similarly in the four cases. However, it is higher for ${\sigma }_r$; Figure 2d shows the ratios between the maximum stresses obtained with FEA in lines 1 and 2 and those obtained with Eq. 8. In Figure 2d, it is noticeable that Eq. 8 offers a better approximation for ${\sigma }_t$, between 0.8973 and 1.0658, at line 1 than for ${\sigma }_r$, between 0.6965 and 1.0133, with respect to the values determined with FEA. For the stresses in line 2, Eq. 8 offers a better approximation for ${\sigma }_r$, with values between 0.9958 and 1.0761, than for ${\sigma }_t$, with values between 0.8589 and 1.0436, with respect to the values determined with FEA. Additionally, Figure 3a shows the ratios of the maximum stresses in Line 2 and the maximum stresses at Line 1. It is noted that the relation between the maximum ${\sigma }_r$ at Line 2 and at Line 1 increases from 1.0248, at $t_{o,i}=1$, to 1.4297, at $t_{o,i}=4.5$. Unlike the relation for ${\sigma }_t$, which decreases from 1.0083, at $t_{o,i}=1$, to 0.9958, at $t_{o,i}=4.5$. Furthermore, Figure 3b and Figure 3c show the results for ${\sigma }_r$, at Line 2, and ${\sigma }_t$, at Line 1, from the FEA compared with the ones obtained with Eq. 8, respectively. Figure 3a and Figure 3b highlight that, although the maximum stresses are fine approximations, the MAPE is considerably high, specially in ${\sigma }_t$ for $t_{o,i}\ge 2.5$. All results from Figure 2 and Figure 3 are presented in function of $t_{o,i}$.

3. Discussion

In this study, the solution of the equilibrium equation for a variable-thickness flywheel, Eq. 3, is based on the assumption that the flywheel thickness varies smoothly, so that $dt/dr\approx 0$, and is given by Eq. 4. This assumption simplifies Eq. 3 and yields Eq. 5. With the assumption that $dt/dr\approx 0$, the solution to Eq. 5 is an easy task. The solution to Eq. 5 is basically the solution for a flywheel of constant thickness plus the term $-BZr/15$, which accounts for the thickness variation. In fact, when the thickness is constant, $B=0$ and $A=t_i$, leading to Eq. 1 and 2.

The validation analysis exhibits a high correlation between the solution, given by Eq. 8, and the results obtained from FEA, for $1\leqq t_{o,i}\leqq 4.5$, Figure 3b. However, P should not be considered as the only way to quantify how good an approximation Eq. 8 is, but rather consider that both solutions, FEA and Eq. 8, behave very much alike within $1\leqq t_{o,i}\leqq 4.5$.

Furthermore, Figure 3c shows that the MAPE between the results, although exhibiting very high values of P, begins to increase as the thickness relation increases. Additionally, in Figure 3d is noted that Eq. 8 offers an excellent approximation for both maximum values of ${\sigma }_r$ and ${\sigma }_t$ within $1\leqq t_{o,i}\leqq 4.5$. Although the maximum ${\sigma }_r$ varies significantly from Line 2 to Line 1, that is, along the thickness of the flywheel, Eq. 8 is an excellent approximation for the maximum ${\sigma }_r$ on Line 2, the line at which the designer is more interested in determining ${\sigma }_r$, since that is where the maximum ${\sigma }_r$ is. Alternatively, Eq. 8 offers a better approximation of maximum ${\sigma }_t$ at Line 1 than at Line 2. Unlike ${\sigma }_r$, the variation of ${\sigma }_t$ along the thickness of the flywheel is too small and could be neglected. Nonetheless, it is important to keep in mind that it does varies and it is higher at Line 1 than at Line 2. However, it can be concluded that Eq. 8 is a very good approximation for maximum ${\sigma }_r$ and ${\sigma }_t$ within $1\leqq t_{o,i}\leqq 4.5$.

Figure 3. (a) Relation of maximum stress at Line 2 and Line 1. (b) ${\sigma }_r$ computed with FEA and Eq. 8. (c) ${\sigma }_t$ computed with FEA and Eq. 8.

Figure 3. (a) Relation of maximum stress at Line 2 and Line 1. (b) ${\sigma }_r$ computed with FEA and Eq. 8. (c) ${\sigma }_t$ computed with FEA and Eq. 8.

Although this solution was developed for t given by Eq. 4, the solution method could be applied to higher-order polynomials of t, as a consequence of Lagrange’s mean value theorem. That is, the solution for maximum stresses would be valid within a 10% variation from FEM as long as $dt/dr\leqq 4.5$.

Equation 8 was observed to approximate ${\sigma }_r$ better than ${\sigma }_t$, Figure 2 and Figure 3. That is attributed to the fact that ${\sigma }_t$ is obtained with the derivative of F, unlike ${\sigma }_r$, which is obtained directly from F. Since F is an approximation, its derivative will lose accuracy and so will ${\sigma }_t$.

Furthermore, although the maximum ${\sigma }_r$ and ${\sigma }_t$ are well approximated, ${\sigma }_t$ exhibits a quantitatively considerable MAPE, higher than 10%, for $t_{o,i}\ge 2.1$, and higher than 15% for $t_{o,i}\ge 3.25$. This is a very important information to bear in mind when using Eq. 8 to design a variable-thickness flywheel. Maximum stresses will be within a fine approximation for high values of $t_{o,i}$, but the higher the thickness relation, the higher a design factor, proportional to the MAPE, should be included in the calculation of ${\sigma }_t$ to avoid under-sizing the flywheel. Finally, it is important to highlight the fact that a high value of $t_{o,i}$ is unusual and the flywheel could be divided into a convenient number of elements and treat every element independently, as recommended by Timoshenko [11].

4. Conclusions

In this paper, the equilibrium equation for a variable-thickness flywheel is approximated analytically, assuming that the flywheel’s thickness varies smoothly, as smoothly as to take $dt/dr\approx 0$. Another assumption was that the thickness varies linearly, that is, given by a first-degree polynomial. However, the solution method could also be applied to polynomials of a higher degree t, as long as $dt/dr$ is negligible or as long as $dt/dr\leqq 4.5$.

The solution shows an accurate approximation of the maximum radial and tangential stress, with respect to FEA, for a wide range of $t_{o,i}$.

Developing this approximation to obtain the stresses in a variable-thickness flywheel gives the mechanical designer an easier, faster, and reliable way to calculate the stresses in the flywheel, without recurring to FEA. Additionally, providing an equation that explicitly relates the stresses in the flywheel and its geometry would allow an easier way to optimize its size, cross section, and topology.