Simulation of SLD Impingement on Wind Turbine Blade Airfoil

Zhu Chengxiang，Zhu Chunling，Fu Bin，Zhao Huanyu.Key Laboratory of Fundamental Science for National Defense-Advanced Design Technology of Flight Vehicle，Nanjing University of Aeronautics and Astronautics，Nanjing 006，P.R.China; .Science and Technology on Space Physics Laboratory，Beijing 00076，P.R.China(Received 8 September 05; revised 0 December 05; accepted 5 January 06)



Simulation of SLD Impingement on Wind Turbine Blade Airfoil

Zhu Chengxiang1*，Zhu Chunling1，Fu Bin2，Zhao Huanyu1
1.Key Laboratory of Fundamental Science for National Defense-Advanced Design Technology of Flight Vehicle，Nanjing University of Aeronautics and Astronautics，Nanjing 210016，P.R.China; 2.Science and Technology on Space Physics Laboratory，Beijing 100076，P.R.China
(Received 18 September 2015; revised 20 December 2015; accepted 5 January 2016)

Abstract:Accurate prediction of droplet impingement on wind turbine blade is one of the most important premises of antiicing and de-icing system design.In a super-cooled large droplets (SLD) conditions，droplet no longer maintains a sphere shape，and it may deform，break up，and splash when moving or impinging on the surface.Semi-empirical models of droplet dynamic behaviours are embedded into the Eulerian droplet model to improve the accuracy of the numerical simulation of droplet impingement limits and local collection efficiency.Eulerian droplet model (Model 1) for small droplets and improved Eulerian droplet model (Model 2) for large droplet are both validated by comparing to the wind tunnel experiment results.Using the proposed methods，droplet impingement limitation and local collection efficiency on the S809 airfoil are calculated in various conditions.A detailed derivation of Model 1 and Model 2 is presented along with a comparison of numerical trajectories，drag coefficient and collection efficiency distributions.The results show that droplet dynamic behaviours，including splashing，break-up and deforming，must be considered to accurately simulate the impingement behaviour in the SLD conditions.And with the increase of the droplet diameters，the effects of the droplet dynamic behaviors on the impingement characteristics are more obvious.

Key words:wind turbine icing; super large droplet; droplet deformation; splash

CLC number: TN925Document code: AArticle ID: 1005-1120(2016)01-0112-09

* Corresponding author，E-mail address: cxzhu@ nuaa.edu.cn.

How to cite this article: Zhu Chengxiang，Zhu Chunling，Fu Bin，et al.Simulation of SLD impingement on wind turbine blade airfoil［J］.Trans.Nanjing Univ.Aero.Astro.，2016，33(1) : 112-120.

http: / /dx.doi.org/10.16356/j.1005-1120.2016.01.112

0　Introduction

Wind turbines are devices harnessing natural wind power to generate electricity.They work in various regions including cold areas，offshore and mountains［1］.Unavoidably，meteorological conditions，containing not only small super-cooled droplets whose diameters are from 10 μm to 50 μmbut also super-cooled large droplets (SLD) whose diameters are larger than 50 μm，are encountered.Then icing on the wind turbine blades takes place.The ice changes the aerodynamic geometry of the wind turbine blades and the load on blades，which affects the aerodynamic forces，leading to significant performance degradation，excessive vibration associated with mass imbalance，and reduction of torque and power output［2-3］.Many anti-icing and de-icing strategies are developed，such as active heating of blades，passive hydrophobic coating on blades，and so on［4］.Wind turbine can operate reliably in icing conditions.In order to design and optimize the ice protection systems，and assess the safety risks arising from wind turbine icing，the procedure of icing，the icing limitations，and heat loads should be actually predicted first.

Numerical simulation methods of ice accretion have been studied for more than 25 years，most of which are used for aircraft icing.Many research institutions have developed their own numerical software，e.g.，LEWICE［5］of NASA，FENSAP-ICE［6］of Numerical Technologies International (NTI )，etc.There aremany differences between aircraft icing and wind turbine icing (for example，the aircraft icing is in general sudden，fatal，and short，while the wind turbine icing is gradual，not fatal but unavoidable and long［2］，with wind turbine blades rotating) the main simulation procedure of ice accretion is similar.Generally，prediction of wind turbine icing also includes three parts.The first part is the simulation of air flow，the second part is the computation of super-cooled droplet trajectories and impingement characteristics，and the last part is calculation of heat and mass transfer on blades.Each part has significant effects on the predict accuracy of ice accretion.Actually，modern computational fluid dynamics (CFD) methods are used to simulate the ice accretion procedure based on the hypothesis that the super cooled droplets are small enough to sustainably maintain sphere.Therefore，if the droplets are small enough，there will be a small error according to the assumption.But when the droplets are larger，for example，with diameters more than 50 μm，they may deform，break up and splash under the action of viscosity force and surface tension，as shown in Fig.1.In these cases，the spherical assumption is no longer valid.

Fig.1　SLD's dynamic behaviors under icing conditions

Lane used Weber number to decide whether the droplet was broken up.Wolfe and Andersen［7］established a droplet splashing model based on the drag.Luliano［8］defined the minus mass source term，which is induced into the SLD splashing model based on the Eulerian droplet model.NASA［9］presented an empirical model of SLD physics on the basis of experimental results.

By referring to the previous research findings on the SLD models，a semi-empirical method is proposed to improve the Eulerin droplet model，so that the droplet motions around the air foil would be simulated and the accuracy of the impingement characteristics would be improved in SLD conditions.SLD's trajectories and local collection efficiency are calculated under different conditions.The results reveal that the changes of the drag coefficient and splashing mass are the reasons for the SLD effects on the impingement characteristics.

1　Mathematical Model

1.1Governing equations

Based on the theory of continuity and momentum conservation，the particular Eulerian formulations of the droplet impingement process in conservative form are presented as

ux，uy，uzare droplet velocities in the directions of x，y and z.L is the length of airfoil chord，and α the volume fraction defined as the ratio of droplet volume in control volumes Vwand VΩ

Here Vnis the droplet velocity in the direction perpendicular to the surface，which can be calculated by the following equation

K in Eq.(1) is the momentum exchange coefficient of air and droplet，which is defined as

where μais the dynamic viscosity of air，ρwthe water density，dwthe droplet diameter.f is the drag function，which can be calculated by

where CDis the drag coefficient and Rewthe droplet Reynolds number of relative velocity.

1.2Droplet deformation

For droplets whose diameters are smaller than 40 μm，it is reasonable to make a hypothesis that the droplets maintain in rigid balls with no deformations throughout the moving process.One of the most important reasons is that forces acted by air on the droplet are small enough to be in balance with surface tension.

In these cases，drag coefficient CD，sphereis calculated in Eq.(6)

In some meteorological conditions，droplets diameters are larger than 40 μm.Air exerts forces on these droplets，and the forces may be too large to break the balance of spherical droplets.Then the droplets may deform to irregular shapes，which in turn affect the air force and the droplet trajectories.

In these cases，drag coefficient CDis calculated in Eq.(7)

where the modified factor b is used to consider the effect of deformation on drag coefficient.It can be obtained using the following equation.

where Weber number We is defined as the ratio of inertia force and surface tension，which can be described as［10］

where ρais the air density，v the air velocity，u the droplet velocity，and σ the surface tension coefficient.

1.3Droplet breakup

Generally speaking，if the Weber number is larger than the critical Weber number［11］，the droplet may break up into many small droplets called the secondary droplets.And the secondary particle size dsecis given by the following equation.

where dw0is the droplet diameter before it is broken up，dsecthe secondary droplet diameter，and μwthe dynamic viscosity of droplet.

1.4Droplet splash

Weber number of droplet impingement Weimpand splash parameter KMare used to determine whether the droplets splash when they collide with the wall［12］.

where unis the droplet velocity in the direction perpendicular to the surface.

According to the experimental data in Ref.［13］，droplets splash when KMexceeds 57.7，and the splash water mass mscan be determined by the following equation［14］

where mimpis the mass of impingement droplet and ω the liquid water content.

2　Model Validation

Actually，there are two mathematical models for droplet impingement prediction in this paper.One model，called the Model 1 is used to calculate droplet trajectories and impingement characteristic under small-droplet-size conditions，which assumes that droplets always keep a sphere shape.Model 1 is established by calculating the governing equations introduced above.The other，called the Model 2 is used to calculate droplet trajectories and impingement characteristic under SLD conditions，which considers the droplet de-formation，break-up and splash.Actually，Model 2 is an improved model on the basis of Model 1.It is obtained by adding the aforementioned deformation model，breaking model and splash model into Model 1to consider the effects of SLD dynamic behaviors on impingement characteristic.

To validate the two numerical simulation methods used in this paper，comparisons between predicted results and IRT experiment data in small droplet condition and SLD condition are shown in Figs.2 (a)，(b)，where horizontal axis is S/c，here，S is the surface distance，c the chord length，and vertical axis parameter β the local collection efficiency.It can be seen that under these two conditions，both the impingement limits and the local collection efficiency predicted in this paper agree well with the experiment results.The mathematical models and numerical simulation methods presented in this paper are reliable.

Fig.2　Comparisons of collection efficiency distributions

3　Results and Analysis

Using the methods proposed in this paper，droplet trajectories and impingement characteristic on a typical airfoil S809 of wind turbine blade are predicted in the conditions，as shown in Table 1.

Table 1　Calculation conditions

Figs.3—6 shows the numerical simulation results in Case 3.The We number distribution around the blade airfoil is shown in Fig.3，which is calculated using Model 2.Fig.4 shows four droplet trajectories.SD1 is calculated using Model 1，and SLD1 is obtained with Model 2.SD1 and SLD1 are droplet trajectories releasing from the same position.For SD2，it is calculated using Model 1，and SLD2，it is obtained with Model 2.SD2 and SLD2 are droplet trajectories releasing from the same position.These two groups of calculation results show that，when droplet diameter is up to 200 μm，there are considerable differences between the droplet trajectories with and with out the consideration of break-up，splash and deformation of droplets.The primary reason is that the drag coefficient and drag force are different，as shown in Fig.5.In Fig.5，the horizontal ordinate is the chord-wise positions，and the vertical ordinate is the calculated drag coefficient using different models (SD1 and SD2 are obtained using Model 1，while SLD1 and SLD2 are gotten with Model 2).It can be seen that when the physicalprocesses of droplet breakup，splash and deformation are considered，and the drag coefficient is much larger than that of the sphere droplet，thus leading to the greater deviation for the droplet trajectories from the air flow.

Fig.3　We number distribution around airfoil (MVD =200μm)

Fig.4　Droplet trajectories around airfoil (MVD =200　μm)

Fig.5　Drag coefficient of droplet around airfoil (MVD =200　μm)

Fig.6　Comparison of local droplet collection efficiency on airfoil (MVD =200　μm)

Comparison of local droplet collection efficiency β using Model 1 (line SD) and Model 2 (line SLD) in Case 3 is shown in Fig.6.The horizontal axis is S/c，where S is the surface distance，c the chord length，and vertical axis is the local collection efficiency.It can be seen that，both impingement limit and mass distribution on the leading surface，which are obtained from simulations in the absence of dynamic actions，exceed the results considering the deformation，splash and break-up of the droplets so much that the effects of super large droplet on the impingement characteristics cannot be ignored.

For other cases，there are similar regulations，and the results are shown in Fig.7.Figs.7 (a)—(g) show the droplet trajectories around the airfoil with different droplet diameters.It can be seen that，with the increasing droplet diameters，the deviations of the two groups of lines (SD1 and SLD1，SD2 and SLD2) are both growing.This attributes to the increasing differences of the drag coefficients，as shown in Fig.8.

Fig.9 indicates that the larger the droplet diameter，the more it affects the impingement characteristics.And when the droplet dynamical behaviors are considered，the mass of the droplet impinging on the surface is reduced.This is not only because of the change of drag force，but also attributes to the splash of droplet.The splash droplet mass in different conditions is shown in Fig.10.

Fig.7　Droplet trajectories around airfoil

Fig.8　Drag coefficient of droplet around airfoil

Fig.9　Comparison of local droplet collection efficiency on airfoil

4　Conclusions

Semi-empirical deformation，breakup and splash models are embed into Eulerian droplet motion model (which is called model1 in this paper) to develop a new droplet motion model (Model 2 in this paper) considering super large droplet dynamic behaviors.Ice accretion on airfoil was simulated by using Model 1 and Model 2.Simulation results are shown to be agreed well with the tunnel experimental results.Themodel simulation reveals that when the SLD impinges the airfoil surface，the effects of the dynamic behaviors on the impingement characteristics cannot be ignored，which is attributes to the change of the drag coefficient and the mass of splash droplets.With the increasing of the droplet diameters，the effects are becoming more and more obvious.So，deformation，breakup and splash of droplets should be considered，especially in the SLD conditions，when ice accretion on the wind turbine blade need to predicted accurately to design efficiency anti-icing or de-icing systems.

Acknowledgements

This work was supported by the National Basic Research Program of China (″973″Program) (No.2014CB046200)，the National Natural Science Foundation of China (Nos.11402114，51506089)，the Fundamental Research Funds for Central Universities，Nanjing University of Aeronautics，and Astronautics Basic Scientific Research Projects for Youth Science (No.NS2014014).

References:

［1］YAN L，KOTARO T，FANG F，et al.A wind tunnel experimental study of icing on wind turbine blade airfoil ［J］.Energy Conversion and Management，2014，85 (9) : 591-595.

［2］RHO M S.Atmospheric icing effects of aerodynamics of

wind turbine blade［C］.Proceedings of the ASME 2013 International Mechanical Engineering Congress and Exposition，IMECE2013-64085.［S.l.］: ASME，2013.［3］YI X，CHEN K，WANG K，et al.Application of CFD technology in wind turbine icing prober design［J］.Transactions of Nanjing University of Aeronautics＆Astronautics，2013(3) : 264-269.

［4］OLIVIER P，ADRIAN I.Anti-icing and de-icing techniques for wind turbines: Critical review［J］.Cold Regions Science and Technology，2011，65(1) :88-96.

［5］WRIGHT W B.Validation results for LEWICE 3.0 ［R］.AIAA 2005-1243，2005.

［6］REID T，BARUZZI G，OZCER I.FENSAP-ICE simulation of icing on wind turbine blades［R］.Part 1: Performance Degradation，AIAA 2013-0750，2013.

［7］WOLFE H E，ANDERSEN W H.Kinetics，mechanism and resultant droplet sizes of the aerodynamic breakup of liquid drops: No.0395-04 (18) SP［R］.California，USA，1964.

［8］LULIANO E，MINGIONE G，PETROSINO F.Eulerian modeling of large droplet physics toward realistic aircraft icing simulation［J］.Journal of Aircraft，2011，48(5) : 1621-1632.

［9］WILLIAM W B，MARK G P.Semi-empirical modeling of SLD physics［R］.AIAA-2004-412，2004.

［10］TAN C S，PAPADAKIS M.Droplet breakup，splashing and re-impingement on an iced air foil［R］.AIAA-2005-5185，2005.

［11］RAIMUND H，WAGDI G H.Eulerian modeling of droplet impingement in the SLD regime of aircraft icing ［R］.AIAA-2006-465，2006.

［12］TRUJILLO M F，MATHEWS W S，LEE C F，et al.Modeling and experiment of impingement and atomization of a liquid spray on a wall［J］.International Journal of Engine Research，2000，1(1) : 87-105.

［13］MUNDO C，TROPEA C，SOMMERFELD M.Numerical and experimental investigation of spray characteristics in the vicinity of a rigid wall［J］.Experimental Thermal and Fluid Science，1997，15(3) : 228-237.

［14］KUOHSING E H，GIAO T V，COLIN S B.Water droplet impingement on simulated glaze，mixed，and rime ice accretions［R］.NASA TM-2007-213961，2007.

Dr.Zhu Chengxiang is a lecturer of the Key Laboratory of Fundamental Science for National Defense-Advanced Design Technology of Flight Vehicle at Nanjing University of Aeronautics and Astronautics (NUAA).Her research interests cover heat and mass transfer，aerodynamic，aircraft icing and ice protection systems，etc.

Dr.Zhu Chunling is a professor and Ph.D.supervisor of the Department of Man-Machine＆Environment Engineering，College of Aerospace Engineering，NUAA.Her research interests include aircraft icing and icing tunnel experiment technology，heat and mass transfer.

Mr.Fu Bin is an engineer of Science and Technology on Space Physics Laboratory.His research interests focus on aircraft icing.

Mr.Zhao Huanyu is a Ph.D.candidate of Department of Man-Machine＆Environment Engineering，College of Aerospace Engineering，NUAA.His research interests lie in numerical simulation of aircraft icing.

(Executive Editor: Zhang Tong)

Fig.10　Mass of splash droplet in Cases 1—8