Preface

Rotating machines are machines that mainly rely on rotating action to complete certain functions of fluid transportation, energy conversion, and power generation. There are many kinds of rotating machines. The typical rotating machines include turbines, compressors, fans, pumps, generators, and motors, which are widely used in electric power, petrochemical, metallurgy, aerospace, and other fields.

The research and application of rotating machines have a long history. In recent decades, energy and environmental issues have become increasingly prominent, and the development and utilization of new energy have been widely valued and developed rapidly. The application of rotating machines in the field of new and renewable energy has become an important development direction. We can often see the beautiful scenery of many wind turbines rotating on the mountains, plains, and sea. This book focuses on the recent research progress and development of rotating machines technology in the energy field. It includes four research reports and one review.

A wind turbine is a device that uses the wind to drive the turbine rotor to rotate to convert the kinetic energy of the wind into mechanical energy, electric energy, and heat energy. Wind energy is one of the fastest-growing and most widely used renewable energies in the world. This book includes three research reports on offshore wind turbines, wind heating by vertical axis wind turbines, and icing impacts on wind turbine blades.

Offshore wind turbines are being increasingly used for generating large-scale wind energy. The application of this technology in the deep sea is an important research field. Factors such as waves, earthquakes, extrema weather, and so on affect the performance of offshore wind turbines. Therefore, their foundation structure is very important. This book compares the structural dynamic response between jacket foundation large-scale offshore wind turbines and monopile ones under wind and seismic loads.

In addition to wind electric power generation, wind heating is also a very important direction of wind energy conversion and utilization. Direct heating of a wind turbine is a form of heating production using wind energy with the advantages of low requirements for wind quality, relatively simple device structure, and high heating efficiency. Vertical axis wind turbines are another important form of a wind turbine. They have characteristics of wind direction independence, novel structure, and easy maintenance, and are suitable for small-scale wind energy utilization. The book includes numerical simulation and wind tunnel tests of heating systems under different working conditions. It also examines the adaptability of the wind turbine and the heating efficiency of the heating device.

Wind turbines sometimes face a variety of extreme weather conditions, such as icing, heavy snow, lightning, sandstorms, and so on. For the wind turbines installed in

cold and humid regions, the icing on the blade surface will often occur, which affects the aerodynamic characteristics and safe operation of the turbine rotor and thus the output performance of electric power of wind farms. This book discusses icing and its impacts on the blade surfaces of wind turbines. It also presents the results of wind tunnel testing of icing on a 2D blade airfoil used for wind turbines. The study tested and analyzed the icing distributions on the blade surface at different attack angles under different ambient temperatures and wind speeds.

Pump turbines are yet another type of rotating machine. A reversible pump-turbine is the core part of a pumped hydro energy storage power station, which can switch between pumping and generating modes based on demand. Signal analysis methods have been widely applied to monitor the condition and fault diagnosis of various kinds of rotating machinery, especially pump turbines. This book reviews signal analysis methods and their applications in reversible pump turbines, including time-domain analysis, time-frequency analysis, mode decomposition, and de-noising. According to this review, signal analysis and processing can help us better understand the operational states of reversible pump turbines. Through mode decomposition of the signal, various components within the signal can be distinguished and signal de-noising can be performed.

A gas turbine is the core component of a power plant. To obtain high efficiency, many new kinds of gas turbines have been proposed. The SOFC-GT hybrid system, which is a combination of solid oxide fuel cells and a gas turbine, is one example. This book presents an exergy analysis of the SOFC-TG system. The exergy balances are established for the different components. According to the results, the exergy efficiencies are improved by increasing the pre-reforming fraction of the air and fuel flows.

We can foresee that rotating machines will play a greater role in more fields in the future. The research on rotating machines will be more in-depth and innovative.

> **Yan Li** Engineering College, Northeast Agricultural University, Harbin, China

**Chapter 1**

Seismic

**Abstract**

seismic

**1. Introduction**

*Kailun Niu and Xinyu Pei*

Comparative Analysis of Transient

Dynamics of Large-Scale Offshore

In this paper, a structural dynamic response comparison between jacket foundation large-scale offshore wind turbines (OWTs) and monopile ones under wind and seismic loads is demonstrated. The interaction between flexible soil and pile foundation is described by Winkler soil-structure interaction (SSI) model. The National Renewable Energy Laboratory (NREL) 5 MW large-scale OWT is studied via the finite element model. The structural transient dynamic response of these two structures under normal operating conditions at rated wind speed and earthquake is calculated. The results show that under the action of seismic and turbulent wind, the jacket has better wind and seismic resistance, and the displacement of the top of the tower is small, which can effectively protect the blades and the nacelle. Compared with monopile, the range of the Mises equivalent stress amplitude of the jacket wind turbine was reduced and the average value was decreased at the time of the seismic. The study also found that the

**Keywords:** OWTs, dynamic response, SSI, jacket foundation, monopile foundation,

The need for renewable resources is becoming more evident as the earth's limited reserves of remaining fossil fuels have essentially halted development and utilization. Among renewable energy sources such as wind energy, solar energy, tidal energy and geothermal energy, wind energy is considered to be the most cost-effective alternative energy and one of the cheapest new power source, showing its great potential to meet demand [1]. The 2021 Global Wind Energy Council (GWEC) released an annual wind report that 2020 was the best year in history for the global wind industry showing

Wind Turbines with Different

Foundation Structure under

*Peilin Wang, Minnan Yue, Chun Li, Yangtian Yan,*

existing jacket design will have a local strain energy surge.

### **Chapter 1**

## Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with Different Foundation Structure under Seismic

*Peilin Wang, Minnan Yue, Chun Li, Yangtian Yan, Kailun Niu and Xinyu Pei*

#### **Abstract**

In this paper, a structural dynamic response comparison between jacket foundation large-scale offshore wind turbines (OWTs) and monopile ones under wind and seismic loads is demonstrated. The interaction between flexible soil and pile foundation is described by Winkler soil-structure interaction (SSI) model. The National Renewable Energy Laboratory (NREL) 5 MW large-scale OWT is studied via the finite element model. The structural transient dynamic response of these two structures under normal operating conditions at rated wind speed and earthquake is calculated. The results show that under the action of seismic and turbulent wind, the jacket has better wind and seismic resistance, and the displacement of the top of the tower is small, which can effectively protect the blades and the nacelle. Compared with monopile, the range of the Mises equivalent stress amplitude of the jacket wind turbine was reduced and the average value was decreased at the time of the seismic. The study also found that the existing jacket design will have a local strain energy surge.

**Keywords:** OWTs, dynamic response, SSI, jacket foundation, monopile foundation, seismic

#### **1. Introduction**

The need for renewable resources is becoming more evident as the earth's limited reserves of remaining fossil fuels have essentially halted development and utilization. Among renewable energy sources such as wind energy, solar energy, tidal energy and geothermal energy, wind energy is considered to be the most cost-effective alternative energy and one of the cheapest new power source, showing its great potential to meet demand [1]. The 2021 Global Wind Energy Council (GWEC) released an annual wind report that 2020 was the best year in history for the global wind industry showing

year-over-year growth of 53%. Installing more than 93 GW wind power in a challenging year with disruption to both the global supply chain and project construction has demonstrated the incredible resilience of the wind industry. According to the GWEC, more than 6 GW of new offshore wind power capacity was installed worldwide in 2020, with China adding more than 3 GW of new offshore wind power capacity, accounting for half of the world's new installed capacity [2]. These estimations reasonably take into account the new improved technologies in wind turbine systems and the offshore wind energy available in considerable quantity. Those sites primarily located at near coastal or offshore areas and experiencing high wind throughout the year, such as the sea off southeast China and the western coast of the United States. However, these areas are with close proximity to the Pacific seismic belts, making them vulnerable to severe and frequent earthquakes [3]. Or rather, OWTs working under such condition generally have to experience dual external loadings primarily stimulated by turbulent wind and seismic excitation. In addition, so as to shrink the cost of energy, the size of OWTs has increased to 9.5 MW recently [4]. Nevertheless, higher wind turbine capacity is accompanied by higher wind turbine towers, which are destined to withstand higher wind pressures and are more vulnerable to seismic loads.

In the past few decades, a large number of studies have been carried out to analyze the dynamic behavior of OWTs under earthquake excitation, in order to improve the stability of wind turbines [5–9], but there is still the issue of oversimplification in the geometry of the model. It is noted that in most of the aforementioned literatures, the rotor and nacelle were either completely ignored or simplified as a lumped mass. In Ref. [10–13], Asareh and Prowell selected FAST as a design basis, developed a seismic module to examine the coupled effect of wind and seismic loads,and employed an improved FAST (also known as NREL Seismic) to investigate the relationship between seismic intensity and structural response. Additionally, in Ref. [14–17], the natural frequency and mode shape of the OWT support structure are calculated using frequency domain, and then the seismic load is calculated according to the seismic response spectrum. However, unsteady aerodynamic loads were often ignored or treated as a rotor thrust, resulting in relatively inaccurate prediction of combined loads acting on blades.

In Ref. [18], the coupled behaviors of seismic and wind loads based on a 2 MW wind turbine was conducted by Witcher. In Ref. [19], Jin et al. also employed the NREL Seismic tool to predict the dynamic responses under multiple hazards associated with earthquake excitation and turbulent wind. In general, the aerodynamic load of large-scale wind turbines increases exponentially with the rotor diameter. Hence, the continuous aerodynamic effects on a larger turbine are determined to be unneglectable in the comparative study of operational and parked states. In Ref. [20], the dynamic response of a 5 MW OWT under seismic loading was further discussed. The finite element model of the SSI was considered for the multi-degree of freedom of fixed foundation. Response spectrum analysis method and transient dynamics analysis method were adopted to verify that the SSI cannot be ignored. In addition, some scholars assumed that there is rigid or linear elastic contact between seabed soil and pile foundation [21–25]. However, due to the severe liquefaction of the seabed soil, the soil porosity and water content change continuously with the pile depth. Furthermore, mechanical dynamic responses on the surface exhibit nonlinear changes in both the horizontal and vertical directions of soil, which can raise the damping of the OWT support structure and lessen the natural frequency. In this case, the relative displacement between soil and structure intensifies the liquefaction of soil during earthquake, *Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

and the soil reaction shows highly nonlinear characteristics related to the structural deformation difference [26–33]. To sum up, an accurate nonlinear SSI model is needed to analyze the dynamic characteristic and response of the OWT support structures under mixed loads.

To achieve the same end, Yang et al. [34] proposed a numerical analysis framework coupled with FAST to obtain more precise responses. However, this method is only applicable to OWTs with monopile support structure, and the results are relatively macroscopic. Yet, it is challenging to examine the local response characteristics of OWTs with different support structures under joint loads.

To tackle these challenges, it is necessary to propose a proper and practical control mechanism of large-scale OWTs under persistent earthquake-and-wind-induced excitations, so as to decrease the volume, prolong the working life, and improve the reliability of the support structures. In the paper, instead of 2 MW OWT, the NREL 5 MW OWT is referred to as the research object, utilizing Winkler soil-structure interaction (SSI) model to describe the interaction between flexible soil and pile foundation. Finite element analysis software ANSYS is used to simulate the dynamic response of OWT foundation under different soils and different seismic intensities, and the influence of seismic load on tower-top displacement and tower-base bending moment is studied. Moreover, the local response characteristics of OWTs with different support structures under seismic load are strictly analyzed. Several principles of nonlinear dynamic response of OWTs structures experiencing sorts of earthquakes are identified. These conclusions will be beneficial in the phase of guaranteeing the optimum performance of OWTs foundations in any operating conditions and in any seismic area.

The remainder of the paper is presented as follows. In Section I, numerical modeling and physical properties are explained. Section II describes wind and seismic loads. In Section III, the finite element method is outlined. Reliability is validated in Section IV. The comparation between the numerical results of two foundations are presented in Section V. Finally, Section VI summarizes the work done and concludes the paper.

#### **2. Numerical simulations**

This paper takes NREL 5 MW OWT as the research object. The support structures of offshore wind power are mainly monopile and jacket, which are usually installed in shallow sea areas. The jacket is mainly composed of grid truss, low steel consumption, convenient transportation and assembly, and good anti-wind and wave performance, suitable for deep water [35]. Since this paper mainly studies the difference between monopile and jacket under the combined action of earthquake and turbulent wind, and wave load has little influence on the jacket, the influence of wave on them is ignored.

#### **2.1 Structural model**

Design parameters were detailed in a report published by the National Wind Technology Center (NWTC), proposed by Jonkman et al. [36]. Besides, the main parameters of NREL 5 MW OWT are shown in **Table 1**.

The wall thickness of the tower decreases linearly with height; the outer diameter of the tower base is 6 m and the wall thickness is 27 mm; the diameter of the tower top connecting the engine room connecting flange is 3.87 m and the wall thickness is 19 mm. A kind of A709 circular section high-strength steel with a density of

#### **Table 1.**

*Main parameters of wind turbine.*

7850 kg/m3 , Young's modulus of 210 GPa, Poisson's ratio of 0.3 and material yield limit of 380 MPa is chosen as the material of tower. So as to take into account the influence of structural paint, flanges, bolts and welding point masses, the calculated density has been modified to 8500 kg/m3 [37]. The structure model and parameters are shown in **Figure 1**.

Although more accurate structural dynamics results can be obtained by using FEM, the FEM method requires high computational resources, so almost all studies ignore or simplify the wind wheel and engine room as the bulk mass. Using a similar simplification, Lavassas et al. [38] studied the structural characteristics of 1.0 MW wind turbine tower under earthquake action, in which the seismic load referred to the relevant Greek seismic codes and guidelines. Hacıefendioğlu [39] used a similar simplified model to study the seismic dynamic characteristics of a 3.0 MW offshore wind turbine under shutdown condition, considering SSI effect. Ma et al. [20] calculated the structural response and load distribution characteristics of the concrete tower of 5.0 MW wind turbine under six different earthquake actions by ABAQUS using a simplified model without considering wind turbines at all.

#### **2.2 SSI modeling**

In this paper, the current Wenkel SSI model of the American Petroleum Institute (API) is adopted to consider the nonlinearity of the lateral and longitudinal soil stiffness [40]. Meanwhile, the SSI model is represented by a nonlinear spring, as shown in **Figure 2**. The calculated monopile and jacket *p-y* curves are illustrated in **Figure 3**.

The burial depth of two structures is 25 m and the burial part is divided into 5 sections, starting from 1.25 m below the seabed, with an interval of 5 m at each calculation step. The nonlinear stiffness of the soil *K*n is based on piecewise calculation, the outside diameter of the jacket pile and the monopile is 2 m and 6 m respectively.

The ultimate bearing capacity of sand varies in depths. The ultimate bearing capacity of shallow and deep sand below the seabed is denoted as Eq. (1) and Eq. (2). *Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

**Figure 1.**

*Structure model and parameters. (a) Wind turbine support structure model diagram, (b) jacket model parameters and (c) jacket finite element model and grid.*

$$p\_s = (\mathbf{C}\_1 \mathbf{H} + \mathbf{C}\_2 \mathbf{D}) \eta H \tag{1}$$

$$p\_{\rm d} = \mathcal{C}\_{\rm 3} D \eta H \tag{2}$$

where *p*s and *p*d are the ultimate bearing capacity of shallow and deep sand respectively; *γ* is the effective gravity of sand; *H* is the depth of sand below the seabed; Coefficient *C*1, *C*2 and *C*3 is determined by the internal friction angle *φ*, displayed in **Table 2**; *D* is the outside diameter of the pile foundation.

**Figure 2.** *Schematic diagram of soil structure interaction.*

When the pile foundation is laterally displaced, the sand reacts to the pile foundation due to the deformation of the pile. The nonlinear relation between the force and the lateral displacement of the sand is denoted as follows:

$$\begin{cases} p = Ap\_s \tanh\left[\frac{kH}{Ap\_s}\mathcal{y}\right] \\ p = Ap\_d \tanh\left[\frac{kH}{Ap\_d}\mathcal{y}\right] \end{cases} \tag{3}$$

*Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

**Table 2.**

*Soil parameters.*

where *A* is the empirical adjustment coefficient with regard to different loads; *A* is 0.9 with cyclic load and *A* ¼ ð Þ 3*:*0 � 0*:*8*H=D* . When with static load, *k* is the ground reaction coefficient; *p* is the sand reaction force and *y* is the lateral displacement.

#### **3. Loading calculation**

#### **3.1 Turbulent wind load**

According to the IEC Kaimal turbulence wind spectrum model defined by the classic IEC61400–1 [41], the three component calculations of the wind (*K=u*,*v*,*w*)are determined by:

$$S(f) = \frac{4\sigma^2 L / \overline{\mathfrak{u}}\_{\text{hub}}}{(\mathbf{1} + \mathfrak{G}\mathbf{L}/\overline{\mathfrak{u}}\_{\text{hub}})^{5/3}}\tag{4}$$

where *f* is the frequency; *u*hub is the mean wind speed at the hub height; *σ* is the standard deviation of the wind speed and *L* is the integral scale parameter of each velocity component. Wind field calculations are arranged at 11 grid points in the y and z directions.

The turbulent component *V t*ð Þ is calculated by applying an Inverse Fast Fourier Transfer (IFFT) to the Kaimal turbulent spectrum described by Power Law wind profile and Logarithmic wind profile as follows:

$$\begin{cases} V(\mathbf{z}) = V(\mathbf{z\_{hub}}) (\mathbf{z}/\mathbf{z\_{hub}})^{0.3} \\ V(\mathbf{y}) = V(\mathbf{y\_{hub}}) \left( \ln \frac{\mathbf{y}}{\mathbf{z\_0}} - \boldsymbol{\nu} \right) / \left( \ln \frac{\mathbf{z\_{hub}}}{\mathbf{z\_0}} - \boldsymbol{\nu} \right) \end{cases} \tag{5}$$

where *V z*ð Þ and *V y*ð Þ represent the wind velocity in the lateral and vertical directions; *V z*ð Þ hub is the mean velocity at the hub height *z*hub; The value of *V z*ð Þ hub is selected as 12.0 m/s equal to the rated wind speed; *z*0 is land surface roughness with a value of 0.021; *ψ* is vertical stability of dimensionless function.

TurbSim [42] developed by NREL is adopted to simulate the full-field turbulent wind. The generated wind field is presented in **Figure 4**. The time-varying wind speed has a peak value of over 20 m/s at the hub and an average magnitude value of 11.4 m/s as expected. The variation of wind speed is irregular in time domain and non-uniform

**Figure 4.** *Turbulent wind at hub height of the OWTs.*

**Figure 5.** *Time history accelerations of ground motions.*

in spatial distribution, which indicates the turbulent characteristics of the generated wind field.

#### **3.2 Seismic load**

The earthquake record was selected from the Pacific Earthquake Engineering Research Center (PEER) NGA database [43]. The earthquake event occurred in Morgan Hill, USA, with a magnitude of 6.19 and a peak acceleration of 2.44 m/s�2 . The time history of seismic acceleration is shown in **Figure 5**.

#### **4. Finite element method**

Transient dynamics analysis is used to analyze the dynamic response of a load structure subjected to any time variation. The Newmark time integral method is adopted to solve the whole structural system matrix for all kinds of large nonlinear elastic–plastic deformation structures [44, 45].

The basic equation of structural transient dynamics analysis is:

$$[M]\{\ddot{u}\} + [C]\{\dot{u}\} + [K]\{u\} = \{F(t)\} \tag{6}$$

*Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with…

The accuracy of the calculation result is improved by using the complete Newmark time integration method. The specific steps are as shown in **Figure 6**.

#### **5. Reliability compliance**

Like any other numerical models, finite element models should be validated before the formal application. Only in this way can we ensure the accuracy of material models, element formulations and mathematical calculations. In this study, the calculation results are compared with those of NREL wind turbine aerodyn-servodyn-hydrodynelasdyn open source software FAST. The dynamic response of tower-top displacement is shown in **Figure 7**.

**Figure 7.** *Time-domain response of tower top displacement under seismic.*

The result from finite element calculation is basically consistent with the result of FAST in the turbulent wind load before 20 s. Both of the curves manifest severe fluctuations, and the trend is highly consistent, with only minor differences in amplitude. The average tower-top displacement is 1.295 m via FAST, but it is 1.234 m using finite element merely with an error value of 4.7%. Therefore, the finite element model is validated.

#### **6. Results and discussion**

#### **6.1 Modal analysis**

In the modal analysis module, Block Lanczos solver performance output is used to extract mode shapes and natural frequencies. The Block Lanczos method uses the sparse matrix solver, overriding any specified solver via the Eqslv command, and can be used for large symmetric eigenvalue problems**.** The 1st and 2nd modal natural vibration frequencies of the two structures involved in SSI are shown in **Table 3**. The 1st and 2nd fore-aft tower vibration mode shapes are illustrated in **Figure 8**.

The vibration characteristics of the OWTs are greatly changed due to the different structures. The 1st natural frequency of the jacket structure is higher, and the 2nd natural frequency is lower, compared with those of the monopile structure. Also, the largest displacement location of the 1st mode at the nacelle is different from that of the 2nd mode.

#### **6.2 Structural dynamic response analysis**

To investigate the effect of wind-seismic interaction on the dynamic responses of the OWTs, a nonlinear time-history analysis is more appropriate compared to the

**Table 3.** *Wind turbine natural frequency.*

*Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

**Figure 8.** *Fore-aft mode shape of tower.*

other analysis schemes. In this study, the dynamic responses of two structures are analyzed under two conditions, with and without seismic load.

Numerical simulations in four diverse cases are carried out:

1. Jacket OWTs under wind and seismic load;

2.Monopile OWTs under wind and seismic load;

3. Jacket OWTs under wind load merely;

4.Monopile OWTs under wind load merely.

Under such circumstance, the simulation duration is 100 s and the time step is 0.002 s. The turbulent wind speed is 11.4 m/s and the earthquake is assumed to occur at 20th s with a magnitude of 6.19. Tower-top displacements in x directions (fore-aft) and y (side-to-side) are presented in **Figure 9**.

When comparing the structure behaviors in different directions, the tower-top displacement in x direction is obviously larger than that in y direction, this is because the x direction suffers from wind load. The displacement curve varies significantly in earthquake, and it is obvious larger in y direction than that in x direction, due to the damping of wind load.

When comparing the different structure behaviors, the dynamic response of the jacket during initial period of the earthquake was smaller than that of monopole, with a 32.3% decrease in peak.

After the earthquake, energy is quickly dissipated in the turbulent wind, the tower-top displacement curve of jacket drops sharply in x direction, with a 39.2% decrease in peak compared with monopile. Most of the energy can be dissipated through the vibration of the jacket tower, which shows that the jacket can be used in aseismic design.

#### **Figure 9.**

*Tower top displacement response in time domain. (a) Earthquake and turbulent wind combined action conditions and (b) turbulent wind action.*

However, the stress of the supporting structure of OWTs will increase sharply when the earthquake occurs, and it is necessary to conduct further research on the response characteristics of shear stress. The average equivalent stress curve of the structures in time-domain is shown in **Figure 10**.

It is shown that the average equivalent stress manifests low-frequency fluctuations as wind speed vary, with a range of 0.92 64.27 MPa for monopile and 0.66 53.31 MPa for jacket. The average value of equivalent stress of jacket decreased *Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

**Figure 10.**

*Average equivalent stress curve in time domain. (a) Earthquake and turbulent wind combined action conditions and (b) turbulent wind action.*

by 21.6% compared with monopile. It is indicated that the jacket structure has a higher and more reliable ultimate bearing capacity of support structures subjected to accumulated fatigue damage.

When seismic load is added, the shear curve shows high frequency fluctuations, accompanied by severe vibration of the tower. The equivalent stress amplitude varies from 1.91 to 82.05 MPa for monopile, and from 5.13 to 59.03 MPa for jacket. The average value of equivalent stress of jacket decreased by 48.7% compared with monopile. Therefore, the jacket structure can reduce the collapse possibility in that increasing rapid stress presumably exceed the tower yield limit when the earthquake hits the coast.

The principle of virtual work states that a virtual (very small) change of the internal strain energy must be offset by an identical change in external work due to the applied loads, and then the virtual strain energy value related to the nodal displacements can be got. The maximum strain energy position of the structure is shown in **Figures 11** and **12** is the cloud picture of strain energy peak moment.

The strain energy fluctuation has been proven to be associated with the change of wind speed, with an average of 228.7 J for monopile and 365.5 J for jacket. Both the two structures receive the same energy from the blades of OWT. However, due to the large displacement response of the monopile, part of the energy is dissipated. The stain energy of stable jacket is 59% higher than that of monopile.

When the earthquake occurred, the local strain energy surged in both structures. The maximum strain energy appeared at 20 s 30 s, the average value of jacket is 1510.8 J, which is 313% higher than that before the earthquake, and the monopile is 708.9 J, which is 94.0% higher than that before the earthquake. Compared with monopile, strain energy of the jacket is larger and more concentrated, which is prone to suffer local deformation.

**Figure 11.**

*Maximum strain energy position curve in time domain. (a) Earthquake and turbulent wind combined action conditions and (b) turbulent wind action.*

**Figure 12.** *Strain energy of OWT foundation structure.*

*Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

#### **7. Conclusion**

In this study, the NREL 5 MW OWT was chosen as the study object and an earthquake with a magnitude of 6.19 is selected as input ground motions. Transient dynamics of large-scale OWTs under combined loads with different foundation structures is analyzed. The structural dynamic response of jacket foundation and that of monopile foundation was carefully compared. In view of the calculations and discussions described, the four key conclusions are as follows:

#### **Acknowledgements**

The authors would like to acknowledge the support of National Natural Science Foundation of China (Grand No. 51976131, 52006148 and 52106262), and this research also supported by the Shanghai "Action Plan for Scientific and Technological" (Grant No.<PHONE_NUMBER>0).

#### **Author details**

Peilin Wang1 , Minnan Yue1 \*, Chun Li1,2, Yangtian Yan1 , Kailun Niu1 and Xinyu Pei1

1 School of Energy and Power Engineering, University of Shanghai for Science and Technology, Shanghai, People's Republic of China

2 Shanghai Key Laboratory of Multiphase Flow and Heat Transfer in Power Engineering, Shanghai, People's Republic of China

\*Address all correspondence to<EMAIL_ADDRESS>
© 2022 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

*Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

#### **References**

[1] Li Y, Sun C, Yu J, et al. Scaling method of the rotating blade of a wind turbine for a rime ice wind tunnel test. Energies. 2019;**12**:626-627

[2] GWEC. Global Wind Report 2021. Brussels, Belgium: GWEC; 2021

[3] Yang Y, Ye K, Li C, et al. Dynamic behavior of wind turbines influenced by aerodynamic damping and earthquake intensity. Wind Energy. 2018;**21**(5): 303-319

[4] Ding Q, Li C, Cheng S, et al. Study on TMD control on stability improvement of barge-supported floating offshore wind turbine based on the multi-island genetic algorithm. China Ocean Engineering. 2019;**33**(3):1-13

[5] Bazeos N, Hatzigeorgiou GD, Hondros ID, Karamaneas H, Karabalis DL, Beskos DE. Static, seismic and stability analyses of a prototype wind turbine steel tower. Engineering Structures. 2002;**24**(8): 1015-1025

[6] Hacıefendioğlu K. Stochastic seismic response analysis of offshore wind turbine including fluid-structure-soil interaction. The Structural Design of Tall and Special Buildings. 2012;**21**(12): 867-878

[7] Zhao X, Maisser P. Seismic response analysis of wind turbine towers including soil-structure interaction. Proceedings of the Institution of Mechanical Engineers, Part K: Journal of Multi-body Dynamics. 2006;**220**(1): 53-61

[8] Sapountzakis EJ, Dikaros IC, Kampitsis AE, Koroneou AD. Nonlinear response of wind turbines under wind and seismic excitations with

soil–structure interaction. Journal of Computational and Nonlinear Dynamics. 2015;**10**(4):041007

[9] Kaynia AM. Seismic considerations in design of offshore wind turbines. Soil Dynamics and Earthquake Engineering. 2019;**124**:399-407

[10] Asareh MA, Prowell I, Volz J, Schonberg W. A computational platform for considering the effects of aerodynamic and seismic load combination for utility scale horizontal axis wind turbines. Earthquake Engineering and Engineering Vibration. 2016;**15**(1):91-102

[11] Asareh MA, Prowell I. Seismic Loading for FAST. (No. NREL/SR-5000-53872). Golden, CO, USA: National Renewable Energy Laboratory; 2011

[12] Yang Y, Bashir M, Li C, et al. Analysis of seismic behaviour of an offshore wind turbine with a flexible foundation. Ocean Engineering. 2019;**178**:215-228

[13] Yang Y, Li C, Bashir M, et al. Investigation on the sensitivity of flexible foundation models of an offshore wind turbine under earthquake loadings. Engineering Structures. 2019; **183**:756-769

[14] Ma H, Zhang D. Seismic response of a prestressed concrete wind turbine tower. International Journal of Civil Engineering. 2016;**14**(8):561-571

[15] Patil A, Jung S, Kwon OS. Structural performance of a parked wind turbine tower subjected to strong ground motions. Engineering Structures. 2016; **120**:92-102

[16] Smith V, Mahmoud H. Multihazard assessment of wind turbine towers under simultaneous application of wind,

operation, and seismic loads. Journal of Performance of Constructed Facilities. 2016;**30**(6):04016043

[17] Tesser RK, Pilla LL, Dupros F, et al. Improving the Performance of Seismic Wave Simulations With Dynamic Load Balancing. Proceedings of the 2014 22nd Euromicro International Conference on Parallel, Distributed, and Network-Based Processing. Los Alamitos, CA: IEEE Computer Society; 2014

[18] Asareh MA, Schonberg W, Volz J. Effects of seismic and aerodynamic load interaction on structural dynamic response of multi-megawatt utility scale horizontal axis wind turbines. Renewable Energy. 2016;**86**:49-58

[19] Jin X, Liu H, Ju W. Wind turbine seismic load analysis based on numerical calculation. Strojniski Vestnik/Journal of Mechanical Engineering. 2014;**60**(10): 638-648

[20] Ma H, Dongdong Z. Seismic response of a prestressed concrete wind turbine tower. International Journal of Civil Engineering. 2016;**14**:561-571

[21] Gao Z, Saha N, Moan T, Amdahl J. Dynamic Analysis of Offshore Fixed Wind Turbines Under Wind and Wave Loads Using Alternative Computer Codes. Crete, Greece: Proceedings of the 3rd EAWE conference, TORQUE 2010: the science of making torque from wind; 2010

[22] Voormeeren SN, van der Valk PLC, Nortier BP, et al. Accurate and efficient modeling of complex offshore wind turbine support structures using augmented superelements. Wind Energy. 2014;**17**:35-54

[23] Banerjee A, Chakraborty T, Matsagar V. Stochastic dynamic analysis of an offshore wind turbine considering

soil-structure interaction. Advances in Structural Engineering. 2014;**10**:673-687

[24] Zhang Y, Liao C, Chen J, et al. Numerical analysis of interaction between seabed and mono-pile subjected to dynamic wave loadings considering the pile rocking effect. Ocean Engineering. 2018;**155**:173-188

[25] Ye J, Wang G. Seismic dynamics of offshore breakwater on liquefiable seabed foundation. Soil Dynamics and Earthquake Engineering. 2015;**76**:86-99

[26] Ibsen LB, Liingaard M. Prototype bucket foundation for wind turbines natural frequency estimation. DCE Technical report No. 9. Aalborg, Denmark: Department of Civil Engineering, Aalborg University; 2006

[27] Anastasopoulos I, Theofilou M. Hybrid foundation for offshore wind turbines: Environmental and seismic loading [J]. Soil Dynamics & Earthquake Engineering. 2016;**80**(1–2):192-209

[28] Galvín P, Romero A, Solís M, et al. Dynamic characterisation of wind turbine towers account for a monopile foundation and different soil conditions. Structure and Infrastructure Engineering. 2016;**13**(7):942-954

[29] Hacıefendioglu K. Stochastic seismic response analysis of offshore wind turbine including fluid-structure-soil interaction [J]. Structural Design of Tall and Special Buildings. 2012;**21**:867-878

[30] Kaynia AM. Seismic considerations in design of offshore wind turbines. Soil Dynamics and Earthquake Engineering. 2018;**124**:399-407

[31] Zhongsheng L, Yang Y, Chun L, et al. Analysis of dynamic seismic response characteristics in timefrequency domain of wind turbine

*Comparative Analysis of Transient Dynamics of Large-Scale Offshore Wind Turbines with… DOI: http://dx.doi.org/10.5772/intechopen.101730*

considering SSI [J]. Journal of Dynamic Engineering. 2018;**38**(7):587-593

[32] Yang Y, Chun L, Quanyong Y. Research on dynamic response of a 5 MW wind turbine tower under seismic conditions. Journal of Chinese Society of Power Engineering. 2017;**11**:83-89

[33] Shen-Haw J, Huang Y-C. Analyses of offshore wind turbine structures with soil-structure interaction under earthquakes. Ocean Engineering. 2019; **187**:106190

[34] Yang Y, Chun L, Weipao M, et al. Dynamic response of wind turbine structure under combined action of turbulent wind field and earthquake excitation [J]. Journal of Vibration and Shock. 2016;**34**(21):136-143

[35] Xiaoni Wu YH, Li Y, et al. Foundations of offshore wind turbines: A review. Renewable and Sustainable Energy Reviews. 2019;**104**:379-393

[36] Jonkman JM, Butterfield S, Musial W, Scott G. Definition of a 5 MW reference wind turbine for offshore system development in: Report No NREL/TP-500-38060. Golden: National Renewable Energy Laboratory; 2009

[37] Agbayani NA. A technical overview of ASCE/AWEARP 2011: Recommended practice for compliance of large landbased wind turbine support structures. Structures Congress. 2014;**2014**: 1759-1770

[38] Lavassas I, Nikolaidis G, Zervas P, et al. Analysis and design of the prototype of a steel 1-MW wind turbine tower. Engineering Structures. 2003; **25**(8):1097-1106

[39] Hacıefendioğlu K. Stochastic seismic response analysis of offshore wind turbine including fluid-structure-soil

interaction. The Structural Design of Tall and Special Buildings. 2012;**21**(12): 867-878

[40] Recommended Practice for Planning, Designing and Constructing Fixed Offshore Platforms-Working Stress Design API. Washington: American Petroleum Institute; 2002

[41] IEC 61400–1. "Wind turbines-Part 1: Design requirements." 3rd ed. Geneva, Switzerland: International Electrotechnical Commission; 2005

[42] Jonkman BJ. TurbSim User's Guide: Version 1.50. (No. NREL/TP-500-46198). Golden, CO, United States: National Renewable Energy Laboratory (NREL); 2009

[43] Pacific Earthquake Engineering Research (PEER) ground motion database. Available from: https://peer. berkeley.edu/

[44] Clough RW. Thoughts about the origin of the finite element method. Computers and Structures. 2001;**79**: 2029-2030

[45] Clough RW. Early history of the finite element method from the view point of a pioneer. International Journal for Numerical Methods in Engineering. 2004;**60**:283-287

#### **Chapter 2**

## Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind Turbine

*Tieliu Jiang, Shengwen Wang, Lidong Zhang and Zhongbin Zhang*

#### **Abstract**

Direct heating of wind turbine is a form of heating production using wind energy. Because of its low requirements for wind quality, relatively simple device structure and high heating efficiency, the wind turbine-driven heating devices can take the place of traditional fossil energy for winter heating and achieve the purpose of reducing carbon emissions to a certain extent. Through the heater experimental platform constructed at Northeast Electric Power University, it is concluded that the liquid stirring heater can operate at 7 m/s, however; the efficiency is extremely low in the low-temperature environment. The permanent magnet eddy current heater must operate at the wind speed above 13 m/s, but it also can operate normally under the low-temperature condition. Combining the advantages and disadvantages of the two experimental devices and setting the vertical axis wind generator as the original motor, the heating efficiency of the two types of heating devices are analyzed under different working conditions, and then the adaptability of the wind turbine and the heating efficiency of the heating device are also studied.

**Keywords:** aerodynamic characteristics, vertical axis wind turbine, wind energy utilization coefficient, wind energy heating system, thermal production efficiency

#### **1. Introduction**

In addition to achieving the energy conversion, wind energy heating system can also get the thermal energy needed by users through the "wind energy-mechanical energy-thermal energy" route [1]. The kinetic energy of the natural wind is captured and converted to mechanical energy by the wind turbine firstly, and then the heater converts the mechanical energy to the desired thermal energy. Compared with the first type of energy conversion, the second form is called the wind energy direct heating system and it saves the power generation equipment (as shown in **Figure 1**)

**Figure 1.** *Schematic diagram of the "wind energy-mechanical energy-thermal energy" conversion pathway.*

and reduces the number of energy conversion time. The system will further reduce the initial investment cost and significantly improve the energy utilization coefficient.

The vertical axis wind turbine has attracted more and more attention due to its simple structure, low cost, and no yaw system required. The vertical axis wind turbine is divided into lift-type and drag-type vertical axis wind turbines, between which the lift-type one has a higher wind energy utilization coefficient under the high blade tip speed ratio, and thus the wind turbine has higher power. Lift-type vertical axis wind turbines are usually designed with two-or three-bladed wind turbines, and the threebladed turbine has a lower shaft torque ripple and better self-starting characteristics compared to the two-bladed turbine [2].

The basic principle of the permanent magnet eddy current heater (**Figure 2**) is that when the rotor of the eddy current heater starts to rotate, the stator heating element is in the changing magnetic field. Since the stator heating element is generally a solid metal structure, many free loops can be formed inside it. Under the action of changing magnetic field, the magnetic flux of each circuit will change, which will produce an induced current. The impedance value formed in the stator heating body is small, so the circuit current will be large, so as to achieve the effect of heating.

**Figure 2.** *Model of the permanent magnet eddy current heater.*

#### *Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

For the study of wind energy heating, the authors in reference [3] set up an outdoor mixing heating experimental platform powered by natural wind and calculated the heating efficiency by the recorded wind speed, rotational speed, and temperature of working fluid. The authors in reference [4] compared the heating effects of the flat blade and the cylindrical blade and found that the heating effect of the flat blade is much better. The authors in reference [5] established a mathematical model to match the torque and the power, which provides a theoretical basis for the design of the stirred wind heating device. The authors in reference [6] optimized the structure of the wind turbine by using the Fluent software according to the relevant theory of wind turbine. The authors in reference [7] used Computational Fluid Dynamics (CFD) method to analyze the thermal efficiency of wind energy heating, and verified the feasibility of using CFD to analyze the mixing heating device.

The basic principle of the liquid stirring heater (**Figure 3**) is that the wind turbine directly drives the agitator to rotate the liquid at high speed and make the liquid heat.

The authors in reference [8] studied the relationship between the torque required in the starting stage of the liquid stirring heater and the stirring impeller radius, angular acceleration. The liquid stirring heater is accompanied by a higher torque when starting.

The authors in reference [9] studied the permanent magnet eddy current heater directly driven by the vertical axis resistance differential wind turbine and analyzed the work of the permanent magnet eddy current heater under a certain wind speed. The authors in reference [10] simulated the permanent magnet eddy current heater model using the finite element method, determined the relevant geometric parameters and material properties of the model and obtained the heater power. The authors in reference [11] set up a permanent magnet eddy current heater experiment device, through which the relevant data were obtained. By using the test device, the temperature changes under different rotating speeds, working times, and different import and export water temperatures were measured, and then the corresponding conversion efficiency was obtained. The authors in reference [12] show that the increase of the thermal energy of the permanent magnet eddy current heater is roughly proportional to the square of the rotational speed increase. The authors in reference [13]

**Figure 3.** *Model of liquid stirring heater.*

**Figure 4.** *Model of vertical axis wind turbine in wind tunnel test.*

pointed out that for the instability and intermittent nature of wind energy, connecting the thermal energy storage device after the permanent magnet eddy current heater can ensure the stable output of thermal energy.

In this study, the operation characteristics of the heater directly driven by a vertical axis wind turbine (**Figure 4**) under different working conditions are studied. The heating efficiencies of the two types of heaters are analyzed, and the matching relationship between the wind turbine and heater is optimized (see **Figures 2**–**4**).

#### **2. Physical model and numerical simulation method**

#### **2.1 Physical model of the wind turbine**

The wind turbine model is a three-blade vertical axis wind turbine adopting the NACA0018 symmetrical airfoil. The string length (*c*) of the blade is 0.25 m, the rotating diameter (*D*) of the wind turbine is 0.55 m. The blade height of the vertical axis wind turbine is 1 m. The geometric model and parameters are shown in **Figure 5** and **Table 1**. The computational domain is mainly divided into two parts, the outside is the static domain and the inner is the internal rotating domain, and sliding grid calculations are applied to rotate the wind turbine in the inner rotating domain. These domains are connected by the interface connection. To eliminate the effect of the inlet blocking on the wind turbine performance and to ensure the continuity of the fluid flow during the simulation, the diameter of the rotating domain inside is set to be 1.5D (see **Figure 6**) [14].

*Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

#### **Figure 5.**

*Geometric model of vertical axis wind turbine.*

#### **Table 1.**

*Dimension parameters of the vertical axis wind turbine.*

**Figure 6.** *Schematic diagram of the 3D Computational domain.*

**Figure 7.** *Force analysis of the element on the blade.*

The aerodynamic characteristics of the vertical axis wind turbine can be expressed as follows:

The following speed relation is satisfied when the wind turbine blades are in each position as Eq. (1):

$$
\overrightarrow{V} = \overrightarrow{U} + \overrightarrow{W} \tag{1}
$$

where, *V* ! —inflow wind speed, m/s; *U* ! —wind turbine blade line speed, m/s; *W* ! relative speed between airflow and blade, m/s.

The force analysis of a certain determined blade is shown in **Figure 7**. The central point of the vertical axis wind turbine is the *O* point, and the height is 2*H*. The *Oxyz* coordinate frame is defined in **Figure 7**, where *Oz* is perpendicular to the vertical rotation axis of the rotor, *Ox* is in the same direction as the wind speed passing through the rotor. And the midpoint of this element in the blade is defined as *M*, string length is *l*.

*W* ! is perpendicular to the blade direction and is expressed as Eq. (2):

$$W\_n = V \sin \theta \cos \delta \tag{2}$$

The other component can be expressed as Eq. (3):

$$W\_t = U + V\cos\theta = r\omega + V\cos\theta\tag{3}$$

So, the force on the blade is expressed as Eq. (4):

$$\boldsymbol{W}\_{u}^{2} = \boldsymbol{W}\_{n}^{2} + \boldsymbol{W}\_{t}^{2} = \boldsymbol{V}^{2}\sin^{2}\theta\cos^{2}\delta + (r\boldsymbol{w} + \boldsymbol{V}\cos\theta)^{2} \tag{4}$$

Based on the above velocity decomposition relationship, the attack angle of blade is as Eq. (5):

$$\tan a = \frac{V \sin \theta \cos \delta}{rw + V \cos \theta} \tag{5}$$

The aerodynamic pressure acting on the blade can be expressed as Eq. (6):

*Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

$$Q = \frac{1}{2}\rho W\_u^2\tag{6}$$

The Lilienthal aerodynamic coefficient at this blade element is expressed as Eq. (7):

$$\begin{cases} \mathbf{C}\_{l} = \mathbf{C}\_{l} \sin a - \mathbf{C}\_{d} \cos a \\ \mathbf{C}\_{n} = \mathbf{C}\_{l} \cos a + \mathbf{C}\_{d} \sin a \end{cases} \tag{7}$$

where *Cl* is the lift coefficient, and *Cd* is the drag coefficient.

Component forces of normal direction and wing string direction of the blade are given as Eq. (8):

$$\begin{cases} dN = \mathcal{C}\_n Q l ds \\ dT = \mathcal{C}\_t Q l ds \end{cases} \tag{8}$$

Decompose the above component forces to the flow wind speed direction, the resultant force received by the rotor in that direction is as Eq. (9):

$$F = \frac{bL}{2\pi} \int\_{-H}^{+H} \int\_{0}^{2\pi} Q(\mathbf{C}\_n \sin\theta \cos\delta - \mathbf{C}\_t \cos\theta) d\theta ds \tag{9}$$

where *b* is the number of blades.

The torque formula provided by the force acting on the blade for the rotor rotation axis is expressed as Eq. (10):

$$dM = dTr = C\_lQLr ds\tag{10}$$

Integrating the above formulas, the torque of the whole rotor is expressed as Eq. (11):

$$M = \frac{bl}{2\pi} \int\_{-H}^{+H} \int\_{0}^{2\pi} C\_l Q r ds d\theta \tag{11}$$

Therefore, the power is expressed as Eq. (12):

$$P = Moo = \frac{bl}{2\pi} \int\_{-H}^{+H} \int\_{0}^{2\pi} \mathbf{C}\_{l} Qroadsd\theta \tag{12}$$

The wind energy utilization coefficient *CP* is expressed as Eq. (13):

$$C\_P = \frac{P}{\frac{1}{2}\rho V^3 A} \tag{13}$$

#### **2.2 Grid division and verification of numerical simulation results**

Grid division is a very important part in numerical simulation. Good grid division can improve the accuracy of the wind turbine performance prediction. In order to ensure the accuracy of the simulation, this study adopts the 3D structured grid for the wind turbine as shown in **Figure 8**.

For the three-blade wind turbine model studied in this paper, the pressure-velocity coupling method and SIMPLE algorithm are used to solve the transient URANS Equation, and the pressure order, momentum term, and turbulence dissipation term

**Figure 8.**

*Mesh of vertical axis wind turbine: (a) top view of mesh in rotating domain; (b) mesh around blade; (c) the entire mesh in the computational domain.*

are all solved by the second-order windward space dissipation algorithm, and the judgment criterion of convergence is set to 105 . The turbulence model used in the simulation is the transition SST turbulence model, which is the original *k-ω*, SST Equations plus two empirical formulas. It is more accurate in the turbulence rotation compared to other simulated wind turbine models. In order to avoid the influences of the grid number, the y+ value, the turbulence model, the algorithm, the time step, and other factors on the simulation results, it needs to verify the simulation results, as shown in **Figure 9**. Compared with the 2D numerical simulation results, it can be

*Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

found that the wind energy utilization coefficient using the 3D numerical calculation results is close to the change trend of the wind tunnel experiment, and the accuracy of the wind energy utilization coefficient prediction is obviously better compared with the 2D numerical simulation. When the wind turbine runs at the optimal tip speed ratio (*λ* = 1.48), the wind energy utilization coefficient obtained from the 3D simulation is 0.297, which is 15% higher than that in the wind tunnel experiment and 19% lower than that in the 2D simulation.

#### **3. Results and discussion**

#### **3.1 Operation characteristics and efficiency of the permanent magnet eddy current heater**

It is known from experiments that when the working fluid temperature of the bulk heat exchange surface increases by 60°C under different speed conditions, the faster the speed of the permanent magnet eddy current heater, the faster the working fluid temperature of the heat exchanger surface increases. However, if the speed is too fast, it will lead to an increase in heat production absorbed by the bulk metal during the heater operation, the uniform heating of circulating working fluid, and a high energy loss. Overall, the running speed of the heater is relatively suitable at 20 rad/s. The starting torque of the permanent magnet eddy current heater is about 6.89Nm, while after the magnetic eddy current heater is used as the wind turbine load, the wind turbine tip speed is relatively low, which will deviate from the optimal tip speed ratio interval, and thus leads to a low wind energy utilization coefficient during operation.

For the permanent magnet eddy current heater, when the wind speed varies from 13 m/s to 17 m/s, the output parameters of the wind turbine are shown in **Table 2**.

The rotating speed and torque curve of wind turbine under different wind speed conditions are shown in **Figures 10** and **11**. Clearly, the rotating speed and torque of the wind turbine both increase linearly with the wind speed increased.

#### **Table 2.**

*Output parameters of wind turbine at 13–17 m/s wind speeds.*

**Figure 10.**

*Curve of the rotating speed change at 13–17 m/s wind speeds.*

**Figure 11.** *Curve of the torque change at 13–17 m/s wind speeds.*

According to the experimental data, the relationship between the heat absorption power *Pa*1 of the permanent magnet eddy current heater during normal operation and the mechanical power *Pw* associated with the input heater is obtained as Eq. (14):

$$P\_{a1} = 0.0008P\_w^2 + 0.3201P\_w + 2.875\tag{14}$$

The relationship between the power of heat absorption *Pa*2 of liquid stirring heater and the mechanical power *Pw* associated with the input heater is obtained as Eq. (15):

$$P\_{a2} = 0.0003P\_w^2 + 0.366P\_w - 10.061\tag{15}$$

*Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

#### **Table 3.**

*The heating efficiency and system efficiency of the permanent magnet eddy current heater.*

Therefore, the efficiency of the heater can be expressed as Eq. (16):

$$
\eta\_{pth} = \frac{P\_a}{P\_m} \tag{16}
$$

where *Pa*—heat absorption power of working fluid, J; *Pm*—input mechanical power of heater, J; *ηpth*—heater efficiency.

The system efficiency is the ratio of the heat obtained by the circulating working fluid to the wind energy swept by the wind turbine. It can directly reflect how much energy the heating system captures from the natural wind. The expression is Eq. (17):

$$
\eta\_{\rm sys} = C\_P \eta\_{\rm path} \tag{17}
$$

According to the above formulas, the heating efficiency and system efficiency of the permanent magnet eddy current heater can be obtained (**Table 3**) by using the output power of the wind turbine.

As seen from **Figures 12** and **13**, with the increase of the test wind speed, the heating efficiency increases at the same time. The heating efficiency and system efficiency are significantly increased when the wind speed is higher than 15.5 m/s. When the wind speed is 17 m/s, the heating efficiency and system efficiency reach the maximum values of 72.21% and 10.83%, respectively. Hence, the permanent magnet eddy current heater has higher efficiency under the condition of higher wind speed and rotating speed.

#### **3.2 Operation characteristics and efficiency of the liquid stirring heater**

According to the experiment results, the temperature rise rate fluctuates in a certain range when the liquid stirring heater rotates at different speeds, but it does not attenuate with the increase of the working fluid temperature. Hence, with the rise of work fluid temperature, the increase of environmental thermal dissipation will not significantly affect the working fluid temperature rise rate. The changes in

**Figure 12.** *Curve of the heating efficiency change at 13–17 m/s wind speeds.*

**Figure 13.** *Curve of the system efficiency change at 13–17 m/s wind speeds.*

working fluid temperature at different heater speeds are in linear function, and the effect of heater speed on heating is significant. The strength and stiffness of mixing blade and the plate under high speed are also need to be considered.

The wind turbine can complete the start-up operation at a low wind speed, and the starting torque is much lower than the permanent magnet eddy current heater at the same power level. Therefore, compared with the permanent magnet eddy current heater, the mixing heater can use the wind energy at a lower speed and improve the wind energy utilization. Due to the small drive torque of the stirring heater, the wind turbine speed is increased, the tip speed ratio is close to the optimal tip speed ratio with a higher wind energy utilization coefficient.

*Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

#### **Table 4.**

*Output parameters of wind turbine at 7–13 m/s wind speeds.*

**Figure 14.** *Curve of the rotating speed change at 7–13 m/s wind speeds.*

When the wind speed varies from 7 m/s to 13 m/s, the output parameters are given in **Table 4**.

The rotating speed and torque curve of the wind turbine in the test wind speed interval are shown in **Figures 14** and **15**. It can be seen that the rotating speed and torque both change linearly in the test wind speed range. Therefore, when the rotating speed of the heater increases, the mixing resistance is borne by the blade and the flow resistance plate will also increase simultaneously. When the heater power is necessary to be further improved, the design of the mixing blade and the damping plate structure should be optimized to increase the mixing resistance and reduce the maximum rotating speed.

Based on Eqs. (15)–(17), the heating efficiency and system efficiency (**Table 5**) of the liquid stirring heater can be obtained according to the output power of the wind turbine.

According to **Table 5**, **Figures 16** and **17**, the heating efficiency increases with the wind speed increased in the test wind speed range, and the maximum efficiency is 48% when the wind speed is 13 m/s. The corresponding wind energy utilization

**Figure 15.** *Curve of the torque change at 7–13 m/s wind speeds.*

#### **Table 5.**

*The heating efficiency and system efficiency of the liquid stirring heater.*

**Figure 16.** *Curve of the heating efficiency change at 7–13 m/s wind speeds.*

*Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

**Figure 17.** *Curve of the system efficiency change at 7–13 m/s wind speeds.*

coefficient is 0.295 with the highest system efficiency. Though the permanent magnet eddy current heater efficiency is high, while the matching characteristics with the wind turbine are poor, leading to low system efficiency. Thus, the good matching of heater and wind turbine will effectively improve the efficiency of the wind energy heating system.

#### **4. Conclusion**

According to the existing experimental results, the numerical simulation method is applied to study the vertical axis wind turbine under different working conditions, the conclusion are as follows:

heater. According to the numerical simulation results, the maximum heating efficiency is up to 72.21%.

3.The heater directly driven by a vertical axis wind turbine system has a certain referred significance for other permanent magnet eddy current heater and liquid stirring heater with vertical axis wind turbine.

### **Author details**

Tieliu Jiang\*, Shengwen Wang, Lidong Zhang and Zhongbin Zhang School of Energy and Power Engineering, Northeast Electric Power University, Jilin, China

\*Address all correspondence to<EMAIL_ADDRESS>
© 2022 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

*Research on Operation Characteristics of Heater Directly Driven by Vertical Axis Wind… DOI: http://dx.doi.org/10.5772/intechopen.102801*

#### **References**

[1] Jianfeng S. Wind energy-thermal energy direct conversion. Cereals and Oils Processing (Electronic Version). 1988;**1**:37-39. Available from: https:// t.cnki.net/kcms/detail?v=8De THzDoIA6bkerXDuqBIgXhgCNDea tEsKvVgKyNfX40D6aYx25vSdkb ssE0Z2uxg64kugubKZyrOEbWori9ubd NtSLb1kWl9Wr2t-u8wdiL0Kps8E3 C3A==&uniplatform=NZKPT [Accessed: January 31, 1988]

[2] Hand B, Kelly G, Cashman A. Aerodynamic design and performance parameters of a lift-type vertical axis wind turbine: A comprehensive review. Renewable and Sustainable Energy Reviews. 2021;**139**. DOI: 10.1016/j. rser.2020.110699 [Accessed: April 1, 2021]

[3] Peng K, Yongguang L. Direct stirring heating experiment powered by the natural wind. Journal of Shanghai University of Electric Power. 2012;**28**(6): 521-524. Available from: https://t.cnki.net/ kcms/detail?v=8DeTHzDoIA7PwdiaE3Hh KPn9DUItEtDrPRzSdSd9PCK0zLkIch JPXjEo7sPs-X4DRb-JXj2kalSK3I951- 9DMbd2cpzCwTLNOPuikETvwl9CLi-UOhceRg==&uniplatform=NZKPT [Accessed: December 15, 2012]

[4] Ting G, Yongguang L, Lihua Z, et al. Experimental study on stirring thermal properties of straight and cylindrical blades. Journal of Shanghai University of Electric Power. 2015;**31**(2):156-160. Available from: https://t.cnki.net/kcms/ detail?v=8DeTHzDoIA6EoFxMlJ8Y39 C3aE\_rTbj0dDDLXp7jjhCOIsjgG\_kLve t1nQhgRXuWJnN-AmXUg7jdphvct60c iADnbHCZ6kNOh66ZiXLDMmSIkH YFCM46f9U6fIh4HGu-&uniplatform= NZKPT [Accessed: April 15, 2012]

[5] Yang L, Yihuai H. Parameter design of stirred wind heating generation device. Acta Energiae Solaris Sinica.

2014;**35**(10):1977-1980. Available from: https://kns.cnki.net/kcms/detail/detail.a spx?FileName=TYLX201410024&DbNa me=CJFQ2014 [Accessed: October 28, 2014]

[6] Lin M. Research on hydraulic wind energy thermal generation system [thesis]. Nanjing: Nanjing University of Science and Technology; 2015. Available from: https://kns.cnki.net/kcms/detail/ detail.aspx?FileName=1016225821.nh& DbName=CMFD2017 [Accessed: December 1, 2015]

[7] Jianzhu Z, Yufeng Z, Guoye W, et al. Study on the thermal efficiency of layered liquid stirring heating system. Acta Energiae Solaris Sinica. 2014;**35**(6): 1034-1039. Available from: https://kns. cnki.net/kcms/detail/detail.aspx?FileNa me=TYLX201406019&DbName= CJFQ2014 [Accessed: June 28, 2014]

[8] Yongguang L, Yunling M. Study on start-up performance of stirring wind thermal device. Acta Energiae Solaris Sinica. 2020;**41**(3):29-33. Available from: https://kns.cnki.net/kcms/detail/detail.a spx?FileName=TYLX202003005&DbNa me=DKFX2020 [Accessed: March 28, 2020]

[9] Fubao W, Maosheng Z, Haipeng T. Analysis of heating generation of vertical axial resistance. Journal of Northwest University (Natural Science Edition). 2013;**43**(4):545-548. Available from: https://kns.cnki.net/kcms/detail/detail.a spx?FileName=XBDZ201304012&DbNa me=CJFQ2013 [Accessed: August 25, 2013]

[10] Nebi O, Fireţeanu V. Finite element electromagnetic 2D model of an eddy current heater with rotating permanent magnets. Annals of the University of

Craiova Romania Series: Electrical Engineering. 2008;**32**(2):62-67. Available from: https://elth.ucv.ro/fisie re/anale/2008/en\_2008\_2.htm [Accessed: October 9, 2008]

[11] Sobor I, Rachier V, Chiciuc A, et al. Small wind energy system with permanent magnet eddy current heater. Telecommunications Energy Conferences Mart Power and Efficiency. 2013;**6**(4):143-150. Available from: http://repository.utm.md/handle/5014/ 14856 [Accessed: November 28, 2013]

[12] Tudorache T, Melcescu L. Outer rotor eddy current heater for wind turbines. Renewable Energy and Environmental Sustainability. 2016;**1**: 25-27. DOI: 10.1051/rees/2016026 [Accessed: June 17, 2016]

[13] Zhaozheng W, Yujun G, Xiaohui Z, et al. Research and development of permanent magnet eddy current heater driven by wind energy. Energy Conservation. 2021;**40**(04):25-28. Available from: https://kns.cnki.net/ kcms/detail/detail.aspx?FileName= JNJN202104008&DbName=DKFX2021 [Accessed: April 25, 2021]

[14] Yan L, Tingting Z, Kotaro T, et al. Numerical simulation on torque characteristics of vertical axis wind turbine with curved plate adding at blade trailing edge. Journal of Northeast Agricultural University. 2017;**48**(3): 72-79. Available from: https://kns.cnki. net/kcms/detail/detail.aspx?FileName= DBDN201703010&DbName=CJFQ2017 [Accessed: March 21, 2017]

#### **Chapter 3**

## Wind Tunnel Test of Icing Distribution on the Leading Edge of a 2D Blade Airfoil for Wind Turbines

*Yan Li, Zhongqiu Mu, Zhiyuan Liu, Wenfeng Guo, Fang Feng and Kotaro Tagawa*

#### **Abstract**

As clean and renewable energy, wind energy has been widely used in the world. The wind turbine is a kind of rotating machinery, which can convert wind energy into mechanical energy and electrical energy. Wind turbines sometimes face a variety of extreme weather conditions, such as icing, heavy snow, lightning, sand storm, and so on, which affect the safety operation of wind turbines. In the present study, icing events on the blade surfaces of wind turbines are focused, and the wind tunnel test of icing was carried out on a 2D blade with NACA0018 airfoil used in wind turbines. In the icing tests, three kinds of ambient temperatures were selected, including −5 °C, −10 °C, and −15 °C, and two kinds of wind speeds were decided, including 5m/s and 10m/s. The icing distributions on the blade surface at the different attack angles were tested and recorded by a high-speed camera for several minutes. The ice accretion and distribution characteristics at the leading edge of the 2D blade airfoil were acquired and analyzed. The findings can provide a reference for the research on the icing mechanism and the de- and anti-icing of wind turbines.

**Keywords:** wind turbine, rotating machinery, blade, icing, wind tunnel test

#### **1. Introduction**

As clean and renewable energy, wind energy has been widely used in the world and will be paid more and more attention [1]. A wind turbine is the most important equipment for the utilization of wind energy, and it is also a kind of rotating machinery, which can convert wind energy into mechanical energy, heat energy, and electrical energy [2]. According to the rotation mode of a wind turbine shaft, there are two key kinds of wind turbines, which are the horizontal axis wind turbine (HAWT) and the vertical axis wind turbine (VAWT). Currently, the HAWTs are used for the large-scale wind turbines mostly. In contrast, the VAWTs are often used for the smallscale wind energy utilization [3]. Thanks to years of researching, the performance of

**Figure 1.** *Icing on the blades of HAWT.*

wind turbines has been greatly improved. The research focus also extends from the improvement of efficiency to how to ensure the safe, stable, and efficient operation of wind turbines. Meanwhile, global climate change intensifies and extreme weather occurs more frequently. Wind turbines are facing more severe weather conditions, such as typhoon, earthquake, lightning, blizzard, icing, etc. [4]. In the cold and moist regions, wind turbines will face the icing problem. When the environment meets a certain condition, supercooled water droplets in the air will accrete on the wind turbine surface and become icing [5]. An icing photo of a wind turbine taken by us in Northeast China in December 2019 is shown in **Figure 1**.

Icing event occurring on the wind turbine blade surface will greatly degrade the efficiency and safety of operation. The research contents of icing mainly include distribution and characteristics of icing on a blade [6–8], aerodynamic performance effect of icing on blade and rotor [9–11], anti-icing and de-icing methods and technologies [12–15], icing detection technologies [16], etc. There are two main icing research methods, numerical simulation [17, 18] and icing wind tunnel test [19, 20]. In the present study, a small-scale icing wind tunnel was used to obtain the distribution characteristics of icing on a model of 2D blade with NACA0018 airfoil. Under different ambient temperatures and wind speeds, experiments on ice accretion on the 2D blade at two kinds of angel of attacks for 4 minutes were carried out. The icing distributions were recorded and analyzed. This research can provide a reference for researching icing mechanism and de- and anti-icing technologies of wind turbines.

*Wind Tunnel Test of Icing Distribution on the Leading Edge of a 2D Blade Airfoil for Wind... DOI: http://dx.doi.org/10.5772/intechopen.104987*

#### **2. Methods**

#### **2.1 Experimental system**

**Figure 2** shows an icing wind tunnel experimental system in Northeast Agricultural University used in this study. It's a low-speed and return-flow type wind tunnel. The refrigeration and spray devices were installed in a normal wind tunnel to provide low temperature and icing environment. The detailed information about this system can be found in the reference [21]. The test section is 250 mm × 250 mm. The range of wind speed is from 1 m/s to 20 m/s. The ambient temperature can be controlled lower to −20°C. By using different sprayers, the icing wind tunnel can supply different kinds of liquid water content (LWC) and medium volume diameter (MVD) of the supercooled water droplet. The LWC is from 0.1 g/m3 ~ 5 g/m3 , and the MVD is from 20 μm ~ 100 μm.

#### **2.2 Test model**

**Figure 3** shows the test model of blade segment with airfoil profile used in wind turbine. To obtain the basic and normal icing characteristics of wind turbine blade, the symmetrical airfoil of NACA0018 was selected because it is often used for wind turbines and basic researches. The blade material was the aluminum alloy, which is usually selected for icing wind tunnel tests because of its good and stable thermal conductivity. According to the size of test section of wind tunnel, the test model has the chord length of 150 mm and the thickness or spanwise length of 20 mm. For this test model, the airfoil profile is consistent along spanwise, and its thickness is small. Therefore, it can be seen as a 2D blade airfoil in the present experimental research. The 3D effect of it can be negligible.

#### **2.3 Test conditions**

Based on the research object, several test parameters were selected for comparison. It included ambient temperatures, wind speeds, and icing times. To check the effect of an angle of attack (AOA) on the blade icing distribution, two kinds of angle

#### **Figure 2.**

*Icing wind tunnel experimental system used in this study. (a) Schematic diagram of the icing wind tunnel experimental system. (b) Photo of the icing wind tunnel.*

#### **Figure 3.**

*Test model of the blade with NACA0018 airfoil.*

#### **Table 1.**

*Experimental conditions of the icing tests in the icing wind tunnel.*

of attacks (AOAs), such as 0°and 10°, were selected. The test conditions are listed in **Table 1**. Furthermore, the LWC in the icing wind tunnel used in the present study was in the range of 0.8 ~ 1.6 g/m3 . The MVD was about 65 μm. In this study, the effects of LWC and MVD on icing of blade airfoil were not researched, which will be carried out deeply in the follow-up study. When the spray system begins to work, a high-speed camera (Phantom v5.1, with the revolution ratio of 1024 × 1024 pixels) begins to record the process of icing on blade surface and take photos at the moments of 2 minutes and the 4 minutes. By using the drawing software, the profile of icing shape can be obtained and analyzed.

#### **3. Results and discussion**

#### **3.1 Icing distributions on the leading edge of blade surface**

**Figures 4** and **5** show the photos of icing distribution at the leading edge of blade under all test conditions. According to the figures, it can be seen that icing occurs on the leading edge of blade for all test conditions. It means that icing will appear when certain environmental conditions are met. However, the characteristics of icing distribution and ice accretion are obviously different under different test conditions.

For the effects of ambient temperature, the most typical characteristic is the change of icing type. When the temperature was −5°C, the type of icing on the blade surface was glaze ice with a transparent color. With the decrease in temperature, the type of rime ice, whose color was white, gradually began to appear. At the temperature of −10°C, there was a little rime ice, which only appeared on the upper surface away from the leading edge. This phenomenon was more obvious at the AOA of 0° than the one at 10°. When the temperature reached to −15°C, most of the ice transformed into the rime ice, while there was very little glaze ice, which only appeared at the front of the leading edge. The glaze ice was just formed by condensation of newly flowing supercooled water droplets in the ultimate period of photography. These

*Wind Tunnel Test of Icing Distribution on the Leading Edge of a 2D Blade Airfoil for Wind... DOI: http://dx.doi.org/10.5772/intechopen.104987*

#### **Figure 4.**

*Photos of icing distribution on blade at the wind speed of 5 m/s.*

**Figure 5.** *Photos of icing distribution on blade at the wind speed of 10 m/s.*

results indicate that the icing type on the leading edge of blade is glaze at high temperatures, such as above −5°C, and transforms into mixed ice and rime ice with the decrease in temperature. Furthermore, another phenomenon should be mentioned. For the AOA of blade at 0°, when the temperature was −5°C, there was a short icicle appearing on the lower surface near the leading edge of the blade at the icing time of 4 minutes. The reason for an explanation of this result is that the heat transfer rate of water droplets decreased with the increase in icing time. At the initial icing stage, the water droplets impinged directly on the blade surface whose temperature was the same as the environmental temperature. In this case, the water droplets froze in a short time because of higher heat conductivity coefficient of aluminum material (237 W/m·K) and low blade surface temperature. That is why the type of rime ice is generated on the blade surface. After that, there was a layer of ice covering the blade surface. The consequent water droplets impinged on the ice layer. The water droplets could not freeze in a short time because of low heat conductivity of ice (2.22 W/m·K). Some super-cooled water droplets attaching on the ice layer surface had not enough time to freeze into ice completely due to the high temperature. These water droplets ran back along the blade surface to the trailing edge, and the icicle was generated. In this case, the type of icing was glaze ice due to high temperature of the icing blade. Under the action of wind speed and gravity, these incompletely frozen water droplets grow and form an ice icicle obliquely downward.

For the effects of wind speed, it was found that the wind speed does not play a key role in the change of icing type. When the temperature was constant, the type of icing was basically the same. When the wind speed increased from 5 m/s to 10 m/s, the icing amount decreases slightly. Meanwhile, the icing shape became smoother, especially at the temperature of −15°C. This may be caused by the constant flow rate of the spray nozzle in unit time. In the present study, the pressure of the pump is constant, which makes the flow rate of the nozzle constant in spite of increasing wind speed. Therefore, the water content passing through a cross-section is constant in unit time. In theory, the icing amounts are the same under two kinds of wind speeds. However, when the wind speed is higher, more water droplets with lightweight are blown away by high-speed wind, and the amount of water droplets impinging on the blade decreases. Additionally, the heat exchange of water droplets with air also accelerates and enhances because of the increase in wind speed. In this case, some small water droplets freeze before impinging on the blade surface, which also results in the decrease in the icing amount.

For the different AOAs, the amount of icing on the lower surface of blade increased because there was an upward angle of attack against wind. In contrast, the amount of icing on the upper surface had a decreasing trend. The key reason is that the windward state of the blade changes from symmetry to asymmetry due to the existence of AOA, which leads to the asymmetry of icing distribution on the leading edge of blade. Therefore, it was concluded that the AOA, or the windward state of blade, determines the location and distribution of icing on blade surface. Furthermore, a deep discussion can be made on how to use these results to the real wind turbine blade. As is known to all, the wind turbine blades always operate in a rotating state. Although the blade model in this research was in a static state, the research findings can be used for the analysis of rotating blades. Based upon the aerodynamics of wind turbine blade, the wind speed and the rotational speed will form a resultant speed at the local element of blade airfoil, which is known as relative wind speed. The AOA in the rotating condition should be decided by the relative wind speed. Therefore, the research findings obtained in the static condition can be used for the analysis of rotating blades.

#### **3.2 Icing profile and icing area**

In order to quantitatively analyze the distribution characteristics of icing, the captured icing photos were processed and the icing contours at the different moments were obtained, which are shown in **Figures 6** and **7**.

For comparatively analyzing with the results at the AOA of 0 degrees, the figures of the icing blade at the AOA of 10° are turned back 10 degrees and located at horizontal level for convenient observation and comparison.

*Wind Tunnel Test of Icing Distribution on the Leading Edge of a 2D Blade Airfoil for Wind... DOI: http://dx.doi.org/10.5772/intechopen.104987*

#### **Figure 6.**

*Ice accretion on blade at the wind speed of 5 m/s.*

#### **Figure 7.**

*Ice accretion on blade at the wind speed of 10 m/s.*

**Figure 8.** *Schematic diagram of icing area.*

Based on these pictures showing the profile of icing, the change of icing shape can be seen clearly. Furthermore, from these figures, the icing area (Ai) is calculated in the present study. The schematic diagram of icing area is shown in **Figure 8**.

**Figure 9** shows the changes in the icing area with the increase in icing time under all test conditions. Additionally, a dimensionless parameter, icing area ratio (*η*s), can also be defined to quantitatively analyze the icing area. It is defined by the ratio of icing area (Ai) and the blade airfoil area (A), which is expressed in Eq. (1).

$$
\eta\_s = \frac{A\_i}{A} \times 100\% \tag{1}
$$

The variations of icing area ratios with icing time, calculated from the data in **Figure 9**, are shown in **Figure 10**.

As shown in **Figure 9**, the icing area increased linearly with the increase in icing time in general. This characteristic can also be found in the result of the icing area ratio in **Figure 10**. Therefore, in this section, it only focuses on the analysis of icing area ratio in **Figure 10**. In the period of 4 minutes icing time, for all conditions in this research, the maximum icing area ratios reached to 6.79% and 8.31% at the AOAs of 0° and 10°, respectively. At the wind speed of 5 m/s, with the decrease in temperature, the growth rate of icing area ratio increased. This phenomenon was obvious, especially at the AOA of 10°. However, for the wind speed of 10 m/s, the effect of temperature on growth rate of icing area ratio became more obvious than that at a wind speed of 5 m/s. The reasons for this result are multifaceted and complex. The high wind speed gives high kinetic energy to the supercooled water droplets, which makes the impact course of water droplets on the blade surface complicated. Also, higher wind speed will affect the characteristics of heat transfer more obviously in comparison with that of the lower wind speed. When the temperature is low, a higher wind speed intensifies the heat transfer of water droplets, which enhances the growth rate of icing area. Deep researches on this issue will be carried out in our further tests.

More discussions were given on the effects of icing time and angle of attack. For the icing time, the icing time length of 4 minutes is selected in this study. Of course, with the increase in icing time, the icing area will increase. In nature, the process of

*Wind Tunnel Test of Icing Distribution on the Leading Edge of a 2D Blade Airfoil for Wind... DOI: http://dx.doi.org/10.5772/intechopen.104987*

**Figure 9.** *Variation of icing area with icing time.*

**Figure 10.** *Variation of icing area ratio with icing time.*

icing on a wind turbine is a long course, which takes time from hours to days. In this research, the main purpose is to explore the initial stage of icing on the leading edge of blade surface. When the ice forms to a certain shape and thickness, the consequent icing process occurs between supercooled water droplets and ice formed on the blade surface, not the blade surface. For the angle of attack, the most important effect is to change the windward area of the blade. Based on the test results of AOA, it is concluded that the icing area will increase along with AOA. For the blade in a static condition, the range of AOA is wide. However, for the rotating blade of a real wind turbine, the relative angle of attack is controlled in a limited range. Therefore, the test results can be a useful reference to the icing research of a real wind turbine.

#### **4. Conclusions**

In this research, the icing tests of a blade model for a wind turbine were carried out in an icing wind tunnel. The distribution characteristics of icing on the leading

#### *Rotating Machines*

edge of the blade were obtained and analyzed. Some key results in the present study are summarized and listed as follows:

Based on the method of icing wind tunnel test in this research, the characteristics of icing on blade surface can be researched and obtained. This study can be a reference to the icing mechanism research and the anti- and de-icing research of wind turbines.

### **Acknowledgements**

This work is supported by the National Natural Science Foundation of China (NSFC) [grant number 51976029]. The authors would like to thanks for the support.

### **Author details**

Yan Li1 \*, Zhongqiu Mu1 , Zhiyuan Liu1 , Wenfeng Guo1 , Fang Feng1 and Kotaro Tagawa2

1 Northeast Agricultural University, Harbin, China

2 Tottori University, Japan

\*Address all correspondence to<EMAIL_ADDRESS>
© 2022 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

*Wind Tunnel Test of Icing Distribution on the Leading Edge of a 2D Blade Airfoil for Wind... DOI: http://dx.doi.org/10.5772/intechopen.104987*

#### **References**

[1] Amjith LR, Bavanish B. A review on biomass and wind as renewable energy for sustainable environment. Chemosphere. 2022;**293**:133579

[2] Li Y. Straight-bladed vertical Axis wind turbines: History, performance, and applications. Rotating Machinery. 2020:87-103

[3] Li Y, Zhao S, Chunming Q, Tonga G, Feng F, Zhao B. Aerodynamic characteristics of straight-bladed vertical Axis wind turbine with a curved-outline wind gathering device. Energy Conversion and Management. 2020;**203**:112249

[4] Dalili N, Edrisy A, Carriveau R. A review of surface engineering issues critical to wind turbine performance. Renewable and Sustainable Energy Reviews. 2009;**13**(2):428-438

[5] Li Y, Wang S, Liu Q, Feng F, Tagawa K. Characteristics of ice accretions on blade of the straight-bladed vertical axis wind turbine rotating at low tip speed ratio. Cold Regions Science and Technology. 2018;**145**:1-13

[6] Wenfeng G, He S, Yan L. Wind tunnel tests of the rime icing characteristics of a straight-bladed vertical axis wind turbine. Renewable Energy. 2021;**179**:116-132

[7] Alessandro Z, Michele DG, Helmut K. Wind energy harnessing of the NREL 5 MW reference wind turbine in icing conditions under different operational strategies. Renewable Energy. 2018;**115**:760-772

[8] Ibrahim GM, Pope K, Muzychka YS. Effects of blade design on ice accretion for horizontal axis wind turbines. Journal of Wind Engineering and Industrial Aerodynamics. 2017;**173**:39-52

[9] Li Y, Shi L, Guo W, Tagawa K, Zhao B. Numerical simulation of icing effect on aerodynamic characteristics of a wind turbine blade. Thermal Science. 2021;**25**(6B):4643-4650

[10] Yan Lie Ce Sun, Yu Jiang, Fang Feng. Scaling method of the rotating blade of a wind turbine for a rime ice wind tunnel test. Energies. 2019;**12**:626-627

[11] Xian Y, Kai-chun W, Hong-lin M, Guo-lin Z. 3-D numerical simulation of droplet collection efficiency in large-scale wind turbine icing. Acta Aerodynamica Sinica. 2013;**31**(06): 745-775 (In Chinese)

[12] Yan L, He S, Wenfeng G. Effect of ultrasonic vibration on the surface adhesive characteristic of iced Aluminum alloy plate. Applied Sciences-Basel. 2022;**12**(5):2357

[13] Yan L, Dong X, Guo W, etc. Ultrasonic micro-vibration deicing of flat plate for wind turbine blades. Journal of Drainage and Irrigation Machinery Engineering. 2022;**40**(2):204-210 (In Chinese)

[14] Li X-j, Guo W-f, Li Y, Zhi X, Feng F. Wind tunnel test of an anti-icing approach by heat pipe for wind turbine blades under the rime ice condition. Thermal Science. 2021;**25**(6B):4485-4493

[15] Zeng J, Song B. Research on experiment and numerical simulation of ultrasonic de-icing for wind turbine blades. Renewable Energy. 2017;**113**:706-712

[16] Xian Y, Yewei G, Guolin Z, et al. Experimental and computational

#### *Rotating Machines*

investigation into ice accretion on airfoil of a transport aircraft. Journal of Aerospace Power. 2011;**26**:808-813 (In Chinese)

[17] Yan L, Long WS, Ce S. Icing distribution of rotating blade of horizontal axis wind turbine based on Quasi-3-D numerical simulation. Thermal Science. 2018;**22**:S681-S691

[18] Yan L, Ce S, Jiang Y, Xian Y, Wenfeng G, Shaolong W, et al. Influence of liquid water content on wind turbine blade icing by numerical simulation. Journal of Drainage and Irrigation Machinery Engineering. 2019;**37**(6):513-520

[19] Li Y, Shen H, Guo W. Simulation and experimental study on the ultrasonic micro-vibration De-icing method for wind turbine blades. Energies. 2021;**14**(24):8246

[20] Guo W, Shen H, Li Y, Feng F, Tagawa K. Wind tunnel tests of the rime icing characteristics of a straight-bladed vertical axis wind turbine. Renewable Energy. 2021;**179**:116-132

[21] Lei S, Fang F, Wenfeng G, Yan L. Research and Development of a smallscale icing wind tunnel test system for blade Airfoil icing characteristics. International Journal of Rotating Machinery. 2021;**2021**:5598859

#### **Chapter 4**

## A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine

*Hao Li, Xianghao Zheng, Yuning Zhang, Jinwei Li and Yuning Zhang*

#### **Abstract**

This chapter summarizes the signal analysis methods and their applications in the field of reversible pump turbine, including time-domain analysis, time-frequency analysis, mode decomposition and de-noising. The time-domain analysis methods include the time-domain statistical parameters, correlation analysis and chaos evaluation. The time-frequency analysis methods include the short-time Fourier transform, Hilbert-Huang transform and VMD-Hilbert transform. Furthermore, the signal mode decomposition and signal de-noising will be introduced together with the related evaluation indexes. The applications of the aforementioned methods are demonstrated based on both the simulated ideal signals and the measured signals from the prototype reversible pump turbines.

**Keywords:** time-domain analysis, time-frequency analysis, mode decomposition, signal de-noising, reversible pump turbine

#### **1. Introduction**

There are a large number of rotating machineries being intensively employed in the modern industry, such as hydro-turbines, pumps, steam turbines, compressors, and generators. They are widely employed in the power generation, heating, refrigeration, and chemical industries. The rotating machineries often involve extreme conditions such as extremely low load, high temperature, high pressure, high speed, and overload. Hence, the malfunctions may occur after the long-term operation with possibly catastrophic consequences. Therefore, it is very important to understand and monitor the operational states of rotating machineries.

With the developments of data analysis technology, signal analysis methods become gradually mature and have been widely applied in the field of condition monitoring and fault diagnosis of various kinds of rotating machineries. The signal analysis could extract useful information from the original signal, and judge the operational states of the equipment based on the obtained information. Signal analysis is of great significance to the rotating machineries both for condition monitoring and fault diagnosis.

A reversible pump turbine is the core part of the pumped hydro energy storage power station. It can switch between the pumping and the generating mode according to the actual demand. When the reversible pump turbine deviates from the design working condition, it is easy to produce serious pressure pulsation and vibration due to the influence of rotor-stator interaction in the vaneless space and the vortex rope in the draft tube, which will affect the normal operation of the pump turbine and the safety of the whole power station [1–3]. Therefore, to monitor the operational states of reversible pump turbine and ensure its safety and efficient operation, it is very necessary to employ the signal analysis methods extensively.

This chapter will review several typical signal analysis methods and introduce their applications in the field of the reversible pump turbine in detail through specific cases. This chapter will be divided into the following four parts. The first part will introduce the time-domain analysis methods of the signals. The second part will introduce the time-frequency analysis methods of the signals. The third part will introduce the signal decomposition and the signal de-noising. The fourth part will demonstrate the applications of the introduced signal analysis methods in the field of the reversible pump turbines with the aid of on-site measured signals.

#### **2. Time-domain analysis methods**

Signals can convey information by expressing the relationship between time and other physical quantities. The time-domain analysis methods of signals refer to a series of processing such as amplification, filtering, statistical feature calculation and correlation analysis in the time domain. Through time-domain analysis methods, the characteristic parameters reflecting the operational state of the mechanical equipment can be extracted from the signal, which could be further employed for the purpose of the evaluation of the operational state and fault diagnosis of the equipment. Time-domain analysis methods could also analyze the auto-correlation (AC) and cross-correlation characteristics of the signal together with the degree of chaos. This section will be divided into three parts including the time domain statistical parameters, the correlation analysis and the chaos evaluation.

#### **2.1 Time-domain statistical parameters**

The classical time-domain statistical parameters include the peak value, peak-topeak value, mean value, mean square value, root mean square, variance and standard deviation. The meanings of these statistical parameters are explained in detail below. Here, the *x t*ð Þ is the signal varying with the time. The *t*1 and *t*2 are the starting and the ending time of the signal in the time domain.

#### **Peak value**

Peak value refers to the maximum amplitude of a signal.

$$\mathbf{x}\_p = \max |\mathbf{x}(t)|\tag{1}$$

#### **Peak-to-peak value**

Peak-to-peak value refers to the difference between the maximum value and the minimum value of a signal within a given period. In the practical application, to eliminate the influences of distortion and noise in the signal, the peak-to-peak value is usually calculated by the confidence interval method with employing the confidence

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

coefficient (e.g., 97 or 95%). Then, the peak-to-peak value could be calculated by using the upper and the lower limits of the obtained signal.

#### **Mean value**

The mean value is the averaged value of the whole signal**.**

$$\overline{\infty} = \frac{1}{t\_2 - t\_1} \int\_{t\_1}^{t\_2} \varkappa(t) dt \tag{2}$$

#### **Mean square value**

Mean square value refers to the mean of the square of the signal, indicating the strength and average power of the signal.

$$\propto\_{ms} = \frac{1}{t\_2 - t\_1} \int\_{t\_1}^{t\_2} \mathbf{x}^2(t)dt \tag{3}$$

#### **Root mean square**

Root mean square refers to the arithmetic square root of the mean square value.

$$\mathcal{X}\_{rms} = \sqrt{\frac{1}{t\_2 - t\_1} \int\_{t\_1}^{t\_2} \varkappa^2(t)dt} \tag{4}$$

#### **Variance**

Variance is the average square value of the difference between the amplitude of the original signal and the mean value. It describes the fluctuation range of the signal and represents the strength of the fluctuating component in the signal.

$$s^2 = \frac{1}{t\_2 - t\_1} \int\_{t\_1}^{t\_2} (\varkappa - \overline{\varkappa})^2 dt \tag{5}$$

#### **Standard deviation**

Standard deviation is the arithmetic square root of variance.

$$s = \sqrt{\frac{1}{t\_2 - t\_1} \int\_{t\_1}^{t\_2} (\varkappa - \overline{\varkappa})^2 dt} \tag{6}$$

#### **2.2 Correlation analysis**

When there are multiple signals, sometimes it is necessary to study the relationship between them. Correlation analysis (including auto-correlation function and cross-correlation analysis) could show the similarity or dependence between signals.

#### **Auto-correlation function**

The auto-correlation (AC) function describes the correlation of a signal between its values at a certain time and after a certain time delay *τ*. It can quantitatively describe the degree of similarity between the time-domain waveform of a signal after shifting *τ* on the time and the original time-domain waveform of a signal.

$$R\_{\mathbf{x}\mathbf{x}}(\tau) = \int\_{-\infty}^{\infty} \mathbf{x}(t)\mathbf{x}(t+\tau)dt = \int\_{-\infty}^{\infty} \mathbf{x}(t-\tau)\mathbf{x}(t)dt\tag{7}$$

*Rxx*ð Þ*τ* is an even function. When *τ* ¼ 0, *Rxx*ð Þ*τ* reaches its maximum value and is equal to the energy of the signal.

#### **Cross-correlation function**

The cross-correlation function describes the degree of similarity of the time-domain waveforms of two signals *x t*ð Þ and *y t*ð Þ at *t* and *t* þ *τ*.

$$R\_{xy}(\tau) = \int\_{-\infty}^{\infty} \varkappa(t)\jmath(t+\tau)dt = \int\_{-\infty}^{\infty} \varkappa(t-\tau)\jmath(t)dt\tag{8}$$

The higher the *Rxy*ð Þ*t* is, the higher the similarity between *x t*ð Þ and *y t*ð Þ þ *τ* is. It should be noted that *Rxy*ð Þ*τ* is not equal to *Ryx*ð Þ*τ* , and the relationship between them is as follows:

$$R\_{\text{xy}}(\tau) = R\_{\text{yx}}(-\tau) \tag{9}$$

#### **Normalization of the correlation function**

In the practical application of signal correlation analysis, to compare the magnitude of the correlation, we often use the normalized form of the correlation function. The normalized form of the auto-correlation function is defined as follows:

$$\rho\_{\rm xx}(\mathbf{r}) = \frac{R\_{\rm xx}(\mathbf{r})}{\sqrt{R\_{\rm xx}(\mathbf{0})}\sqrt{R\_{\rm xx}(\mathbf{0})}} = \frac{R\_{\rm xx}(\mathbf{r})}{R\_{\rm xx}(\mathbf{0})} \tag{10}$$

The normalized form of the cross-correlation function is defined as follows:

$$\rho\_{\rm xy}(\mathbf{r}) = \frac{R\_{\rm xy}(\mathbf{r})}{\sqrt{R\_{\rm xx}(\mathbf{0})}\sqrt{R\_{\rm yy}(\mathbf{0})}} \tag{11}$$

where *ρxx* j j ð Þ*τ* ≤1, *ρxy*ð Þ*τ* � � � � � �≤1. When *ρ τ*ð Þ¼ 1, it means that the correlation between the two signals reaches the maximum; when *ρ τ*ð Þ¼ 0, it means that two signals are not correlated; when *ρ τ*ð Þ¼�1, it means that two signals are linearly negatively correlated.

#### **2.3 Measurement of chaos**

Chaos is an inherent characteristic of a nonlinear dynamic system. When a signal is generated by a nonlinear system, it may show the characteristic of chaos. At this time, the signal waveform is very irregular and is very similar to the characteristic of random noise. Entropy can be adopted to measure the chaos of a nonlinear dynamic system. The more chaotic the system, the higher its entropy. Therefore, we can use entropy to evaluate the degree of chaos in a signal. As a typical method, the permutation entropy for evaluating the degree of chaos in a signal is introduced below [4].

#### **Permutation entropy**

Permutation entropy is an index to measure the complexity of the time series [4]. The more regular the time series is, the smaller the permutation entropy is. The corresponding calculation process is given as follows [4].

Considering a time series *x* ¼ *xi* f g , *i* ¼ 1, 2, … , *N* , it can be reconstructed into the following form [4]:

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

$$\begin{cases} X\_1 = \left\{ \mathbf{x}\_1, \mathbf{x}\_{1+\tau}, \dots, \mathbf{x}\_{1+(em-1)\tau} \right\} \\ \vdots \\ X\_i = \left\{ \mathbf{x}\_i, \mathbf{x}\_{i+\tau}, \dots, \mathbf{x}\_{i+(em-1)\tau} \right\} \\ \vdots \\ X\_{N-(em-1)\tau} = \left\{ \mathbf{x}\_{N-(em-1)\tau}, \mathbf{x}\_{N-(em-2)\tau}, \dots, \mathbf{x}\_N \right\} \end{cases} \tag{12}$$

Here, the "*N*" is the total number of the data points within the whole signal. The "*τ*" refers to the delay time and the "*em*" refers to the embedded dimension.

Now, each *Xi* could be sorted in the ascending order internally [4]:

$$\mathfrak{x}\_{i+\left(\begin{array}{c}j\_1-1\end{array}\right)\mathfrak{r}} \leq \mathfrak{x}\_{i+\left(\begin{array}{c}j\_2-1\end{array}\right)\mathfrak{r}} \leq \dots \leq \mathfrak{x}\_{i+\left(\begin{array}{c}j\_m-1\end{array}\right)\mathfrak{r}} \tag{13}$$

Here, the " *j* 1, *j* 2, … , *j em*" in the subscripts refer to the indexes for different columns. Thus, a sequence of symbols can be obtained for any vector *Xi* [4]:

$$S(h) = \begin{bmatrix} j\_1, j\_2, \dots, j\_{em} \end{bmatrix} \tag{14}$$

There are up to "*em*!" kinds of such index sequences. The occurrence probability of each index sequence is calculated and expressed as *P*1, *P*2, … , *Pr*, respectively. Here, *h* ¼ 1, 2, … ,*r* and *r*≤*em*!. Then, the permutation entropy of the original signal can be expressed as follows [4]:

$$H\_p(em) = -\sum\_{h=1}^r P\_h \ln P\_h \tag{15}$$

When *Ph* ¼ 1*=em*!, the above formula can reach the maximum value ln ð Þ *em*! . Therefore, the normalized form of Eq. (15) can be obtained as follows [4]:

$$H(em) = \frac{H\_P(em)}{\ln\left(em!\right)}\tag{16}$$

The variation range of the normalized permutation entropy is between 0 and 1. In the practical application of calculating the permutation entropy of a signal, the whole data points within the signal can be intercepted by a window with a certain length. And the permutation entropy of the data points in each window is calculated by moving the window one data at a time. In the present chapter, the average permutation entropy of all the windows will be taken as the permutation entropy of the whole signal.

#### **2.4 Examples of a signal time-domain analysis**

In this section, the pressure pulsation signal (the pressure pulsation signal 1) in the vaneless space of a prototype pump turbine under the dimensionless load condition (the ratio of the actual operational load to the rated load) *P\** =52% is taken in the analysis to demonstrate the applications of the time-domain analysis method introduced in Sections 2.1–2.3. The length of the signal is 1 s and the sampling frequency is 2000 Hz. In this section, *p\** refers to the dimensionless pressure pulsation using the

water head (326.01 m). **Figure 1** shows the dimensionless pressure pulsation 1 in the time domain.

#### **Time-domain statistical parameters**

The statistical parameters of the dimensionless pressure pulsation signal 1 in the time domain are shown in the **Table 1** as follows:

#### **Auto-correlation analysis**

**Figures 2** and **3** show the normalized AC functions of the dimensionless pressure pulsation signal 1 and a random noise signal. When *τ* ¼ 0, the AC function of the dimensionless pressure pulsation signal 1 in the **Figure 2** reaches its maximum value, but it gradually decreases to near 0 on both sides of *τ* ¼ 0. For the random noise signal, its AC function reaches the maximum also at *τ* ¼ 0. When *τ* 6¼ 0, due to the

#### **Figure 1.**

*The time-domain diagram of the dimensionless pressure pulsation signal 1.*

#### **Table 1.**

*The time-domain statistical parameters of the dimensionless pressure pulsation signal 1.*

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

**Figure 2.** *The normalized AC function of the dimensionless pressure pulsation signal 1.*

**Figure 3.** *The normalized AC function of the random noise.*

randomness of noise, its AC function decreases rapidly. Therefore, the normalized AC function of the random noise in **Figure 3** is close to 1 at *τ* ¼ 0 and approximately equal to 0 on both sides of *τ* ¼ 0.

#### **Cross-correlation analysis**

In the experiment, another dimensionless pressure pulsation signal (the dimensionless pressure pulsation signal 2) was measured in the vaneless space at the same

**Table 2.**

*Analysis results of permutation entropy of different types of signals.*

**Figure 4.**

*The time-domain diagram of the dimensionless pressure pulsation signal 2.*

operational condition (*P\** =52%) during the different time period. **Figure 4** shows the time-domain diagram of the dimensionless pressure pulsation signal 2.

When the delay time *τ* ¼ 0, the normalized cross-correlation function of the dimensionless pressure pulsation signals 1 and 2 is 0.98, indicating that there is a strong correlation between them. On the contrary, when the delay time *τ* ¼ 0, the normalized cross-correlation function of the dimensionless pressure pulsation signal 1 and random noise signal is 0.0088, indicating that there is almost no correlation between them.

#### **Permutation entropy analysis**

In the present chapter, the embedded dimension, delay time and window length required in the analysis procedure of the permutation entropy are 6, 3 and 128, respectively. The analysis results of the permutation entropy of the dimensionless pressure pulsation signal 1, the dimensionless pressure pulsation signal 2 and the random noise signal are as follows:

It can be seen from the **Table 2** that the analysis results of the permutation entropy of the dimensionless pressure pulsation signals 1 and 2 are smaller than that of the random noise, because the dimensionless pressure pulsation signals contain certain

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

characteristic frequencies with regularity. The permutation entropy of the random noise is large, which indicates its high degree of chaos.

#### **3. Time-frequency analysis methods**

The aim of the time-frequency analysis is to study the time-varying signals, which can reflect the relationship of the variations of the frequency and amplitude with time within the signal. Generally speaking, the typical time-frequency analysis methods include the short-time Fourier transform [5], the Hilbert-Huang transform [6] and the VMD-Hilbert transform with a detailed introduction in this section.

#### **3.1 The traditional Fourier transform**

Before introducing the time-frequency analysis methods, the traditional Fourier transform is:

$$F(\alpha) = \int\_{-\infty}^{\infty} \varkappa(t) e^{-j\alpha t} dt \tag{17}$$

where *x t*ð Þ is the original signal and *e*�*jωt* is the basis function of the Fourier transform. The traditional Fourier transform can accurately reflect the frequency components contained in the signal, but it cannot reflect the time information of these frequency components.

#### **3.2 The time-frequency analysis methods**

#### **Short-time Fourier transform**

$$\text{STFT}\_{\mathfrak{x}}(\tau,\omega) = \int\_{-\infty}^{\infty} \mathfrak{x}(t)\mathfrak{g}(\tau-t)e^{-j\omega t}dt\tag{18}$$

The formula of the short-time Fourier Transform (STFT) is presented above [5]. Compared with the traditional Fourier transform, the window function is added to the formula. STFT needs to select a window function, which can be moved continuously to identify the frequencies at different times. The window function will affect the resolution of the spectrum. The smaller the window size, the higher the time resolution of the spectrum. Normally, the low frequency band requires high frequency resolution and low time resolution. And, the high frequency band requires low frequency resolution and high time resolution. However, after the window function is selected, the resolution of the entire spectrum is fixed and cannot be adjusted. Therefore, STFT is suitable for analyzing the segmented stationary signals or approximately stationary signals.

#### **Hilbert-Huang transform**

Hilbert-Huang transform (HHT) consists of two parts including the empirical mode decomposition (EMD) and the Hilbert spectrum analysis (HSA). These two parts will be introduced in detail below.

The EMD method can decompose the original signal into several intrinsic mode functions (IMFs) and the residual component *r t*ð Þ. The specific steps are given as follow.

Firstly, calculating all the extreme points of the signal *x t*ð Þ, and using the cubic spline function to interpolate and fit the maximum and minimum points respectively to obtain the upper envelope *a t*ð Þ and lower envelopes *b t*ð Þ of the signal;

Secondly, calculating the mean curve *c t*ð Þ of the upper and lower envelopes and the difference *h t*ð Þ between the signal *x t*ð Þ and *c t*ð Þ [6]:

$$c(t) = \frac{a(t) + b(t)}{2} \tag{19}$$

$$h(t) = \mathfrak{x}(t) - \mathfrak{c}(t) \tag{20}$$

Thirdly, the screening stop condition is given as [6]:

$$\sum\_{t=0}^{T} \frac{\left| h\_{\mathbf{g}-1}(t) - h\_{\mathbf{g}}(t) \right|^2}{h\_{\mathbf{g}-1}^2(t)} \le \varepsilon \tag{21}$$

Here, the "*g*" refers to the number of screenings. The "*T*" refers to the signal time length and the "*ε*" refers to a preset value. If *h t*ð Þ does not meet the condition, replacing*x t*ð Þ with *h t*ð Þ and repeating the above steps. And if *h t*ð Þ meets the condition, then *h t*ð Þ will be the first IMF to be extracted.

Fourthly, repeating the above steps to obtain all the remaining IMF until the residual component *r t*ð Þ becomes a constant value sequence with a frequency approximately 0.

Although, the EMD can realize signal decomposition, it has the problems of endpoint effect and mode mixing. Endpoint effect means that analysis results at the signal endpoints will produce errors, which will affect the internal data. Mode mixing means that the EMD method cannot effectively separate different mode components based on the time characteristic scales. These problems can lead to serious performance degradation of EMD methods.

HSA performs a Hilbert transform on each IMF obtained by EMD. The specific steps are given as follows:

Firstly, performing the Hilbert transform on each IMF *hi*ð Þ*t* [6]:

$$H[h\_i(t)] = h\_i(t) \* \frac{1}{\pi t} = \frac{1}{\pi} \int\_{-\infty}^{\infty} \frac{h\_i(\tau)}{t - \tau} d\tau \tag{22}$$

Here, *i* represents the *i-*th IMF component; "\*" represents the convolution operation.

Secondly, using the above formula to construct the analytical signal *di*ð Þ*t* of each IMF *hi*ð Þ*t* [6]:

$$d\_i(\mathbf{t}) = h\_i(\mathbf{t}) + jH[h\_i(\mathbf{t})] = a\_i(\mathbf{t})e^{j\rho\_i(\mathbf{t})} \tag{23}$$

Here [6],

$$a\_i(t) = \sqrt{h\_i^2(t) + H^2[h\_i(t)]}\tag{24}$$

$$\rho\_i(t) = \arctan\left\{\frac{H[h\_i(t)]}{h\_i(t)}\right\} \tag{25}$$

Thirdly, the instantaneous frequency of each IMF can be obtained from the following formula [6]:

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

$$w\_i(t) = \frac{d\rho\_i(t)}{dt} \tag{26}$$

After expressing the results on the time-frequency plane, the Hilbert amplitude spectrum of each IMF *hi*ð Þ*t* can be obtained as follows [6]:

$$H\_i(w, t) = \begin{cases} \mathbf{0}, w \neq w\_i(t) \\ a\_i(t), w = w\_i(t) \end{cases} \tag{27}$$

Fourthly, by presenting the amplitude spectrum of all IMFs in a spectrum, the Hilbert spectrum of the original signal can be obtained as follows [6]:

$$H(w, t) = \text{Re}\sum\_{i=1}^{n} a\_i(t)e^{j\int w\_i(t)dt} \tag{28}$$

The advantage of HHT is that it is not restricted by the linear and stationary characteristics of the signal and can analyze non-linear and non-stationary signals. At the same time, the HHT is self-adaptive through using the EMD to eliminate the choice of basis function. In addition, HHT is not restricted by Heisenberg's uncertainty principle.

#### **VMD-Hilbert transform**

The principle of VMD-Hilbert transform is similar to HHT. It also uses the signal decomposition method to obtain the IMF component of the signal, and then obtains the Hilbert spectrum of the signal through HAS. The difference between these two is that the variational mode decomposition (VMD) [7] is adopted as the signal decomposition method.

The core idea of the VMD method is to construct and solve the variational optimization problem in the frequency domain. The specific steps include the following ones.

Firstly, decomposing the original signal into *K* IMFs with respective central angular frequencies and limited bandwidths [7]:

$$\text{IMF}\_k(t) = A\_k(t)\cos\left[\rho\_k(t)\right] \tag{29}$$

Here, *φk*ð Þ*t* refers to phase of the signal and *Ak*ð Þ*t* refers to instantaneous amplitude.

Secondly, the constraints are added to ensure that the sum of the estimated bandwidths of each IMF is minimum and the sum of *K* IMFs is equal to the original signal [7].

$$\begin{cases} \min\_{\{\text{IMF}\_k\}, \{aw\_k\}} \left\{ \sum\_{k=1}^K \left\| \partial\_t \left[ \left( \delta(t) + \frac{j}{\pi t} \right) \* \text{IMF}\_k(t) \right] e^{-jaw\_k t} \right\|\_2^2 \right\} \\\\ \text{s.t.} \sum\_{k=1}^K \text{IMF}\_k(t) = \mathbf{x}(t) \end{cases} \tag{30}$$

Here, *δ*ð Þ*t* refers to the Dirac function. *ωk* is the central angle frequency of the *k*-th IMF*k*. "\*" refers to the convolution operation.

Thirdly, the Lagrange multiplication operator is introduced to transform the constrained variational problem into an unconstrained variational problem. And an augmented Lagrange expression is obtained as follows [7]:

$$\begin{aligned} L(\{\text{IMF}\_{k}\}, \{\boldsymbol{\nu}\_{k}\}, \boldsymbol{\lambda}) &= \alpha \sum\_{k=1}^{K} \left\| \boldsymbol{\vartheta}\_{t} \left[ \left( \boldsymbol{\delta}(t) + \frac{\boldsymbol{j}}{\pi t} \right) \* \text{IMF}\_{k}(t) \right] \boldsymbol{e}^{-j\boldsymbol{\alpha}\_{k}t} \right\|\_{2}^{2} + \\ \left\| \mathbf{IMF}(t) - \sum\_{k=1}^{K} \text{IMF}\_{k}(t) \right\|\_{2}^{2} &+ \left< \boldsymbol{\lambda}(t), \text{IMF}(t) - \sum\_{k=1}^{K} \text{IMF}\_{k}(t) \right> \end{aligned} \tag{31}$$

Here, *α* refers to the quadratic penalty factor.

Fourthly, the alternating direction multiplier method is used to iteratively solve IMF*k*, *ωk* and *λ* to find the saddle point of the unconstrained variational problem. The update formulas of IMF^ *kn*þ1ð Þ *ω* , *ωn*þ1 *k* and ^*λ n*þ1 ð Þ *ω* obtained by Parseval Fourier equidistant transformation could be expressed as [7]:

$$\text{I\hat{M}F}\_{k\_{n+1}}(o) = \frac{\hat{\mathfrak{x}}(o) - \sum\_{i=1}^{k-1} \text{I\hat{M}F}\_{i}^{n+1}(o) - \sum\_{i=k+1}^{K} \text{I\hat{M}F}\_{i}^{n}(o) + \frac{\hat{\lambda}^{\text{n}}(o)}{2}}{1 + 2a(o - o\_k^n)^2} \tag{32}$$

$$\rho\_{k\_{n+1}} = \frac{\int\_0^\infty \alpha \left| \mathbf{I} \hat{\mathbf{M}} \mathbf{F}\_k(\alpha) \right|^2 d\alpha}{\int\_0^\infty \left| \mathbf{I} \hat{\mathbf{M}} \mathbf{F}\_k(\alpha) \right|^2 d\alpha} \tag{33}$$

$$
\hat{\boldsymbol{\lambda}}^{n+1}(\boldsymbol{\alpha}) = \hat{\boldsymbol{\lambda}}^{n}(\boldsymbol{\alpha}) + \chi \left( \hat{\boldsymbol{\kappa}}(\boldsymbol{\alpha}) - \sum\_{k=1}^{K} \mathbf{I} \hat{\mathbf{M}} \mathbf{F}\_{k}^{n+1}(\boldsymbol{\alpha}) \right) \tag{34}
$$

Here, "^" refers to the Fourier transform operation. The "*γ*" refers to the noise tolerance. And the "*n*" refers to the number of iterations in the calculation process.

When the calculation result meets the given solution accuracy *ε*, the iterative operation is stopped and *K* IMFs are outputted [7].

$$\sum\_{k=1}^{K} \left( \frac{||\mathbf{I}\hat{\mathbf{M}}\mathbf{F}\_k^{u+1} - \mathbf{I}\hat{\mathbf{M}}\mathbf{F}\_k^u||\_2^2}{||\mathbf{I}\hat{\mathbf{M}}\mathbf{F}\_k^u||\_2^2} \right) < \varepsilon \tag{35}$$

Fifthly, the IMF^ *kn*þ1ð Þ *ω* is transformed from the frequency domain to the IMF in the time domain through the operation of the inverse Fourier transform.

Compared with the EMD, the VMD overcomes the problem of mode mixing and has a more solid mathematical theoretical foundation.

In the present chapter, the "*K*" is 7, the "*α*" is 2000, the "*γ*" is 0 with "*ε*" is 10�7 .

#### **3.3 Examples of signal time-frequency analysis**

This section uses the vibration signal at the top cover of a prototype reversible pump turbine during the start-up transient process of generating mode as the analysis object to demonstrate the applications of the time-frequency analysis methods introduced in Sections 3.1 and 3.2. The length of the signal is 2 s and the sampling frequency is 1000 Hz. **Figure 5** shows the time-domain diagram of the vibration signal at the top cover.

**Figure 6** shows the result of the fast Fourier transform (FFT) of the vibration signal. It can be seen from the **Figure 6** that the most obvious frequency component in the vibration signal is 40 Hz with an amplitude of about 0.46 mm. In addition, there is also a relatively obvious peak near 80 Hz with an amplitude of about 0.17 mm.

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

**Figure 5.** *The time-domain diagram of the vibration signal at the top cover.*

However, **Figure 6** cannot provide information about the amplitude variation of different frequency components with time, which is a shortcoming of the traditional Fourier transform.

**Figures 7** and **8** are the STFT results of the vibration signal with the window length of 64 and 256 respectively. From these two figures, not only the amplitudes of different frequency components in the vibration signals but also the variations of the amplitude of each frequency varying with the time can be seen. However, the STFT results are significantly affected by the window size. By comparing **Figures 7** and **8**, it can be concluded that there exist obvious amplitudes at and around 40 Hz in the **Figure 7**, while it almost has obvious amplitudes only at 40 Hz in the **Figure 8**, which is caused by the difference in the window size between them. In the actual application of the STFT, it is difficult to choose the appropriate window size, which is the main problem of this method.

Different with the fixed basis functions of traditional Fourier transform and STFT, HHT and VMD-Hilbert transform can obtain basis functions through adaptive mode decomposition. In addition, HHT and VMD-Hilbert transform do not need to select the window functions and they can achieve high resolutions both in time and frequency domains. **Figures 9** and **10** are the results of the HHT and VMD-Hilbert

**Figure 7.** *The STFT result of the vibration signal with the window length of 64 and the window function of Blackman.*

**Figure 8.**

*The STFT result of the vibration signal with the window length of 256 and the window function of Blackman.*

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

**Figure 9.** *The analysis result of HHT of the vibration signal.*

#### **Figure 10.**

*The analysis result of VMD-Hilbert transform of the vibration signal.*

transform of the signal, respectively. It can be found that the analysis result in the **Figure 9** is unclear because the dominant frequencies of the signal cannot be clearly observed, while the analysis result in the **Figure 10** is very clear. This can be attributed to that the VMD overcomes the mode mixing problem in the EMD adopted in the HHT.

#### **4. Signal mode decomposition and signal de-noising**

Mode refers to the natural vibration characteristics of the mechanical structure, which is composed of natural frequency, damping ratio and mode shape. In practical engineering, the vibration of rotating machinery is often composed of different complex modes. The process of decomposing the vibration signal into sub-signals with different frequency bands on the basis of different modes is called the mode decomposition of the signal. Through the mode decomposition, the modes containing noise

can be found and will be further deleted through the signal de-noising. In our previous research [8], similar methods with those shown in this section have been adopted to analyze the vibration signal from the rotating machinery of nuclear power. However, the specific frequency components between the signals between the rotating machinery of nuclear power and the reversible pump turbine are a little different.

#### **4.1 Mode decomposition of simulated swing signal**

Typical signal analysis methods, such as EMD and VMD, have been fully described in Section 3.2. In this section, a simulated swing signal at the turbine guide bearing of a prototype pump turbine will be used to demonstrate the signal de-noising effects based on these two methods. The simulated rotational speed of the unit is 500 rpm. The characteristic frequencies contained in the simulated signal are 1 time, 2 times and 3 times, the rotational frequency of the impeller with the corresponding amplitudes of 10, 2 and 1 μm, respectively. The signal is dimensionless by dividing the corresponding amplitude of the rotational frequency of the impeller (10 μm). The length of the signal is 2 s and the sampling frequency is 1000 Hz. The white noise with the signal-noise ratio (*SNR*) of 20 dB is added to the simulated swing signal. **Figure 11** shows the time-domain diagrams of the simulated swing signals before and after

**Figure 11.**

*The time-domain diagrams of simulated swing signals. (a) The simulated swing signal without noise. (b) The added noise with the* SNR *of 20 dB. (c) The simulated swing signal with noise added.*

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

adding the white noise. **Figures 12** and **13** are the EMD and VMD analysis results of the simulated swing signal, respectively.

#### **4.2 Results of signal de-noising**

Through mode decomposition, several IMFs are obtained to further constitute the given signal with low noise. The IMFs contain different frequency bands. By removing

**Figure 12.** *The EMD analysis results of the simulated swing signal.*

**Figure 13.** *The VMD analysis results of the simulated swing signal.*

the IMFs dominated by the random noises and retaining the IMFs dominated by the useful signals, the signal de-noising will be proceeded. Generally speaking, the crosscorrelation coefficient between the noise-dominated IMF and the original signal is small, while the cross-correlation coefficient between the IMF component dominated by the useful signal and the original signal is relatively large. Based on this principle, the noise-dominated IMFs can be identified according to the cross-correlation coefficients.

**Table 3** shows the cross-correlation coefficients of the IMFs obtained by EMD and VMD and the simulated swing signal. Here, the threshold of the cross-correlation coefficient is set to 0.12 [8]. And the IMFs with the cross-correlation coefficients larger than this threshold are retained.

Based on the analysis results in **Table 3**, the IMF3 and IMF4 in the EMD analysis results and the IMF1 and IMF2 in the VMD analysis results are selected to reconstruct the de-noised swing signals, respectively. **Figure 14** shows the de-noising results of the simulated swing signal based on EMD and VMD, respectively.

**Figure 14** shows that both the de-noised signals based on EMD and VMD are relatively smooth oscillation curves. Compared with the simulated swing signal without noise, the de-noised signal based on EMD has some abnormal spikes or irregularities as indicated by the red circles in the **Figure 14b**, while the time-domain waveform of the de-noised signal based on VMD is much closer to the simulated swing signal without noise, which indicates that the signal de-noising based on VMD has a better effect than that based on EMD.

#### **4.3 Evaluation indexes of signal de-noising**

The evaluation indexes of signal de-noising can help us quantitatively analyze the de-noising effects. Some evaluation indexes will be introduced below. Here *x t*ð Þ refers to original simulated swing signal without noise; "*i*" refers to the *i*-th data point of signal *x t*ð Þ which contains *M* data points in total; *x*~*i*refers to de-noised swing signal.

**Signal-noise ratio (SNR)**

$$\text{SNR} = \mathbf{10lg}\left(\frac{\sum\_{i=1}^{M} \mathbf{x}\_i^2}{\sum\_{i=1}^{M} (\mathbf{x}\_i - \tilde{\mathbf{x}}\_i)^2}\right) \tag{36}$$

The larger the SNR, the better the de-noising effect.

**Table 3.**

*Cross-correlation coefficients of IMFs based on EMD and VMD and the simulated swing signal.*

*A Review of Signal Analysis Methods and Their Applications in the Reversible Pump Turbine DOI: http://dx.doi.org/10.5772/intechopen.102883*

**Figure 14.**

*The de-noising results of the simulated swing signal based on EMD and VMD. (a) The simulated swing signal without noise. (b) The de-noised signal is based on EMD. (c) The de-noised signal is based on VMD.*

#### **Correlation coefficient (R)**

$$R = \frac{\sum\_{i=1}^{M} \mathbf{x}\_i \tilde{\mathbf{x}}\_i}{\sqrt{\sum\_{i=1}^{M} \mathbf{x}\_i^2 \sum\_{i=1}^{M} \tilde{\mathbf{x}}\_i^2}}\tag{37}$$

The larger the *R*, the better the de-noising effect. **Root-mean-square error (RMSE)**

$$RMSE = \sqrt{\frac{1}{M} \sum\_{i=1}^{M} \left(\mathbf{x}\_i - \mathbf{\tilde{x}}\_i\right)^2} \tag{38}$$

The smaller the RMSE, the better the de-noising effect.

The evaluation indexes of signal de-noising effects based on EMD and VMD are calculated respectively and the results are shown in **Table 4**.