CrossMark

Available online at www.sciencedirect.com

ScienceDirect

Energy Procedía 49 (2014) 1582-1591

SolarPACES 2013

Numerical simulation of wind loads and wind induced dynamic

response of heliostats

C.C. Zanga*, J. M. Christian13, J. K. Yuanb, J. Smentb, A. C. Moyab, C. K. Hobf, and Z.F.

aKey Laboratory of Solar Thermal Energy and Photovoltaic System of Chinese Academy of Sciences; Institute of Electrical Engineering Chinese

Academy of Sciences, 6 Beiertiao Zhongguancun, Beijing 100190, China b Concentrating Solar Technologies Department, Sandia National Laboratories, P.O. Box 5800, Albuquerque, NM87185-1127, USA

Abstract

This paper presents a numerical simulation method for wind loads fluctuations on heliostats and analyzes the dynamic response of the heliostat which is compared with the tested data. Based on the tested wind data, the target spectrum of wind velocity was determined for the numerical simulation of wind velocity. The time history of wind velocity fluctuation was simulated using an Auto-Regressive (AR) mathematical model. Under the simulated wind velocity fluctuation, a modal and transient analysis was performed using the computational model of the heliostat in ANSYS to determine the dynamic response of the heliostat structure under the transient loads. In this paper, the directional deformation in the direction of wind load was obtained under the wind pressure fluctuation. Then, the deformation was converted to acceleration by differentiation using MATLAB. The simulated time history of acceleration was analyzed in the frequency domain and compared with the tested acceleration. For the simulated acceleration, the dominant frequency value of 3.16 Hz had the highest probability of occurrence. For the tested acceleration, the corresponding frequency was 3.18 Hz. These results demonstrate the wind-load simulation method and dynamic response analysis presented in this paper for heliostat performance analysis. In future work, methods for wind-load simulation or calculation should be explored to further improve the simulation accuracy.

©2013The Authors. Published by ElsevierLtd.This is an open access article under the CC BY-NC-ND license (http://creativecommons.Org/licenses/by-nc-nd/3.0/).

Selection and peer review by the scientific conference committee of SolarPACES 2013 under responsibility of PSE AG. Final manuscript published as received without editorial corrections.

* Corresponding author. Tel.: +86-10-82547037; fax: +86-10-62587946.

E-mail address: zangchch@mail.iee.ac.cn t Corresponding author. Tel.: +1-505-844-2384; fax: +1-505-845-3366. E-mail address: ckho@sandia.gov

1876-6102 © 2013 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http:// creativecommons. org/licenses/by-nc-nd/ 3.0/).

Selection and peer review by the scientific conference committee of SolarPACES 2013 under responsibility of PSE AG.

Final manuscript published as received without editorial corrections.

doi:10.1016/j.egypro.2014.03.167

C.C. Zang et al. /Energy Procedia 49 (2014) 1582 - 1591 Keywords: Heliostat; Numerical simulation of wind load; Dynamic response; Dominant frequency

1. Introduction

The major goal of this work is to obtain the dynamic response of heliostats under wind load excitation, which may provide references to optical performance analysis and structural optimization for heliostats. For the solar power tower station, there are always hundreds of heliostats which track the sun and concentrate the sunlight into the receiver fixed on the tower. Once the deformation of the heliostat structure exceeds the allowable value, the beam image will be probably enlarged or move out of the receiver. Wind loading is the main factor that may cause the deformation. Therefore, it is important to research on the wind loading and wind-induced dynamic performance of heliostats. In recent years, some achievements on this subject have been acquired. Bo. [1, 2] et al measured the wind pressure on the surface of the heliostat through wind tunnel test and calculated the wind-induced displacement, equivalent stress and the natural vibration frequency of the heliostat. Huss [3] et al. developed a novel methodology for evaluating the dynamic effect of wind load on Ivanpah LH-2 heliostats. They employed an aero-elastic heliostat model for wind tunnel test to measure the total wind forces and moments which reflect the total fluctuating forces including the background loading as well as the inertial loading. Therefore, the peak dynamic loading could be obtained. Wind tunnel test is an effective way to study wind load on heliostats. Besides that, wind loading test on a real heliostat and numerical simulation of wind load are worth of exploration.

In this paper, a numerical simulation method of wind load and wind-induced dynamic response analysis are performed on the heliostat which is located at Sandia National Laboratories' National Solar Thermal Test Facility (NSTTF), Albuquerque NM, USA. For the aim of low cost and high quality of optical characterization, a lot of studies have been done on the NSTTF heliostat. Ho [4] et al. did dynamic testing and analysis of the heliostat. The objectives of these tests and analysis are to characterize and understand the impacts of dynamic wind load on heliostat optical performance and structural fatigue. Griffith [5] et al. studied the method of structural dynamic measurements for design evaluation and monitoring of heliostat vibrations. They tested vibration modes, strain and displacement under wind loading, which is used to evaluate and improve structural models that predict the motions or deformations of the heliostat due to gravitational and dynamic wind loading. This paper studies the numerical simulation of wind load and analyzes the dynamic response of the heliostat under the wind load excitation. This work aims to provide a method of wind load calculation for heliostat designers.

Nomenclature

a Undetermined parameter

C Spatial correlation matrix

Damping matrix

F (t ) External loads on the structure

G System gain

K Karman's constant

[K] Stiffness matrix

m Length of the sequence

[M] Mass matrix

n Frequency

N (t ) Normal random number, mean value is zero and standard deviation is 1

P Order of AR model

R C/At) The cross-correlation function of two points at time (/At)

S(n) Power Spectral Density (PSD) of the wind velocity fluctuation

S (n) The power spectrum density of wind velocity fluctuation

At Time step

u (t) Wind speed fluctuation

u f. Friction speed

{u} Displacement vector of nodes

{u} Velocity vector of nodes

{ti} Acceleration vector of nodes

U (z) Mean wind speed at the height of z

v Mean wind speed

vf (t) Time history of independent wind speed fluctuation

V (t) Wind velocity vector

w( m) White noise excitation

x( m) Output signal

z Height over the ground

Z0 Roughness length

(pk AR model parameter (A: = 1,2...p)

Standard deviation of the pure random process

2. Approaches

This paper covers four parts of work: 1) wind velocity and direction as well as acceleration measurement; 2) numerical simulation of wind load fluctuation on the heliostat; 3) dynamic response analysis of the heliostat under simulated wind load; 4) simulated and tested results comparison. Based on the measurement, the wind profile near the heliostat was obtained, which was used to simulate the wind load fluctuation on the heliostat. Under the simulated wind load, the dynamic response, such as deformation and acceleration, was analyzed in ANSYS. The simulated acceleration was finally compared with tested acceleration.

2.1. Tests on NSTTF heliostat

To match the simulated wind load with the real wind load on heliostats as well as possible, the wind velocities at different heights and the corresponding wind direction need to be measured. Wind profile could be obtained from the wind velocity testing values and then the wind spectrum target for the numerical simulation of wind load was chosen. Based on the wind direction and the orientation of tested heliostat, the direction of simulated wind load on ANSYS model could be determined. In addition, the acceleration at different locations of the heliostat was tested to evaluate the simulated acceleration.

2.1.1. Description of the test on NSTTF heliostat

The wind velocity and direction was measured by three ultrasonic anemometers mounted on a tripod which was located in the unobstructed flat land west of the instrumented heliostat. These three anemometers were respectively mounted at 20, 25, and 30 feet above the heliostat's base. Based on the wind velocities at different heights, the wind profile around the heliostat was obtained using mathematic method of curve fitting. A power law was used to estimate the wind profile which describes the approaching wind velocity (V) at any height (h). The power function was defined as V=2.77h0 57. The wind profile is shown in Fig.1.

Fig. 1 Wind profile

Elevation Motor

Torque Tube

otaranl-

..Reference DC Accelerometer

DC Tri-Axial Accelerometer

Strain Gauge Data Acquisition

Fig.2 Heliostat sensors layout

Fig.3 Accelerometer number

In order to perform the dynamic response test under wind load excitation, 24 tri-axial DC MEMS accelerometers with sensitivities of 1000 mV/g were mounted via magnets at the yoke and the ends of the trusses which support the facets of the heliostat, as indicated by the red blocks in Fig.2. The 1000 mV/g accelerometers were selected based on previous preliminary modal analysis which found this sensitivity to be optimum. The accelerometer number is shown in Fig.3.

The heliostat positions are defined by the elevation and azimuth angle. Three positions of heliostat were tested for the wind loading and acceleration measurement which covered elevation angle of 45 °and azimuth angle of 180 ° 225 ° and 270 ° respectively. In this paper, the analysis for the position of 180 ° azimuth angle was performed. Due north was defined as the origin for the wind direction and the azimuth angle of the measured heliostat. The azimuth angle of 180 ° means the heliostat faces due south.

2.1.2. Testing data processing

Based on the testing data of wind velocity, the Power Spectral Density (PSD) of the wind velocity fluctuation was obtained and compared with three standard wind velocity spectrums, shown in Fig.4. The three spectrums [6] include modified Kaimal spectrum, Antoniou & Asimakopoulos spectrum and Davenport spectrum which are described respectively in Equation (1), (2) and (3).

nS(n) 100f

nS (n)

nS (n)

(0.44 + 33f )5/

"(0.44 + 5/ )5

(1 + x 2)4/3

KU (z)

ln(—)

X = ■

U(Z) 1200«

It can be seen from Fig.4 that the tested wind velocity spectrum agrees with the standard spectrums. In this paper, Davenport's spectrum is chosen to be the target spectrum of the wind simulation for the following reasons: a) Davenport's spectrum is widely used and proved its reasonableness and reliability. b) The expression of Davenport's spectrum is relatively simple and its assumed condition applies to this paper. c) Fig.4 shows the tested wind spectrum agrees well with the Davenport's spectrum.

2.2. Numerical simulation of dynamic wind loads on heliostat

PSD of tested wind velocity Davenport spectrum Kaimal spectrum Anto & Asi spectrum

10 10 Frequency (Hz)

Fig.4 PSD Comparison

To determine the time history of wind load is the basis for analyzing the dynamic response of heliostat. Generally, the wind load on structures can be obtained using wind tunnel test, field measurement and numerical simulation by computer. Either wind tunnel test or field measurement is relatively not economy and takes long time, while the numerical simulation may overcome that disadvantages. Therefore, the numerical simulation of wind load was performed in this paper.

Wind velocity fluctuation changes irregularly with time. 1 The statistical data of wind velocity demonstrates that the time history of wind velocity fluctuation basically follows stationary Gaussian random processes. Methods for 1 numerical simulation of wind velocity fluctuation are mainly classified into two types: harmony superposition method and linear filtering method, which are also

respectively called spectral representation method and time series method. Auto-Regressive (AR) model is one of the valid methods of linear filter to simulate wind velocity time history and has been widely applied on random vibration and time series analysis due to the advantages of small computational amount and fast computational speed. Any stationary random signal may be considered as the output signal when white noise excites a stable invertible causal system. By using the AR model and linear filter technology, the white noise excitation with zero average value is transferred into stable random signal with specified spectrum characteristics [7, 8]. The AR model is described in Equation (7).

x(m) + ^ a{x(m - i) = Gw(m) (7)

At any moment, the wind velocity vector V (t) is composed of mean wind speed v and fluctuating wind speed u (t) . The mean wind speed may be calculated according to the wind profile near the heliostat.

V (t ) = v + u(t) (8)

The fluctuating wind speed vector u(t )has the character of spatial correlation and is calculated using equation (9).

u (t ) = C *v(t) (9) The vector v(t) is made up of a series of independent fluctuating wind speed time history vt (7) which is calculated using the AR model.

Time History of Wind Velocity

Time History of Wind Velocity

60 Time t(s)

60 Time t(s)

Fig. 5 Simulated time history of wind velocity at facet 5 (Left) and facet 25 (Right)

v (t )=Ew (t -k At )+°nN (t)

(pk and On are calculated by Equation (11)-(13).

Д, (jAt) = ¿Д. [(7 - k) At] cph и = 1,2... P)

R ( jAt) = jSv (n) • cos(2^n • jAt) • dn

= R (0Mv (к At)

(11) (12) (13)

This paper presents the simulation of the time history of wind velocity on NSTTF heliostat using said AR model. The mean wind velocity at different heights was determined by the test wind profile. The wind velocity fluctuation was simulated using the AR model in MATLAB. As seen in Fig.4, Davenport spectrum may be determined as the power spectrum for the simulation. According to the equations (9)-(13), the wind velocity fluctuation time history at the center of 25 facets of NSTTF heliostat was obtained respectively. Since the effect of the spatial correlation on the wind velocity was considered, the time history of wind velocity fluctuation is related to the orientation of heliostats. In this paper, the wind load time history was analyzed at multiple orientations of the heliostat. Fig.5 shows the time history of wind velocity at facet No.5 and No.25 when the heliostat is at 45 degrees in elevation angle and 180 degree in azimuth. The facets number is shown in Fig.6. The simulated sample time step is 0.1 s and total steps number is 1200.

The simulated time history of wind velocity was discretized and converted to power spectral density (PSD) using the mathematical statistics method in MATLAB. To validate the accuracy of the simulation, the simulated PSD of wind velocity was compared with the target PSD of Davenport spectrum. Fig.7 displays the comparison results. It can be seen that the simulated PSD agrees

Fig.6 Facet number of NSTTF heliostat

Facet 5

Fig.7 Power spectral density comparison at facet 5 (Left) and facet 25 (Right) with the target Davenport's spectrum. The wind load simulation is validated.

2.3. Dynamic Analysis of Heliostats

To obtain the dynamic response of the heliostat under wind load excitation, this paper presented the modal and transient analysis using ANSYS. The directional deformation in the wind direction was obtained through the analysis and was converted to acceleration through differentiation so as to be compared with the tested acceleration.

2.3.1 ANSYS model of NSTTF heliostat

The NSTTF heliostat geometry was created in SolidWorks. The geometry represents the NSTTF field heliostat 11-West-14 and was canted to the 180 foot level on the NSTTF tower. The geometry was composed of surfaces for all applicable components and solids for the remaining components. Surfaces allow for FEA model simplification by reducing the overall mesh element count by using shell elements. The geometry was imported into the ANSYS WorkBench mesher. The surface components were meshed with shell elements and the solid components were meshed with tetrahedral elements. 1682850 elements were needed to provide a grid-independent solution. Most of these elements were located in the small geometric solid features which require many elements to accurately solve the problem physics. The mesh was then imported into ANSYS Mechanical to solve the dynamic heliostat problem.

A modal-to-transient linked solution system was used to solve the FEA. The linked system approach first solves a modal analysis providing the heliostat mode shapes and natural frequencies. The solution of this system is linked to the input parameters of a transient analysis. The transient analysis solves for displacements as a result of time varying input forces. The linked system approach allows for quicker computations than a pure transient analysis will provide. Another advantage of the linked system is that for each displacement in time, the accompanying mode shape/frequency for each loading condition can be output.

The modal and transient studies share the same geometry, mesh, and boundary conditions. Three important boundary conditions

Fig.8 NSTTF Heliostat Boundary Conditions, (Left-Zoomed) Fixed torque tube BC, (Middle-Zoomed) Fixed heliostat base BC, (Right-Zoomed) Revolute joint on torque tube BC

„ — Simulation data for facet 23 _ Testing data at 20 feet high

- / \ 1 M / ^ V 1 ' ir

7 / / M 1 i*M

f /¥i /

- fWil .

Simulation data for facet 25 Testing data at 20 feet high

¡vy+^llJ

> 7.5 -

1.6 1.8 Time (s)

1.6 1.8 Time (s)

Fig.9 Simulated and tested wind velocity comparison at facet 23 (Left) and facet 25 (Right)

are applied in this model. The heliostat base is fixed in spaced using a fixed constraint (no translation or rotation). The heliostat torque tube has two different boundary conditions applied. The fixed end of the torque tube was connected to the elevation motor which was assumed to be rigidly fixed (no translation or rotation). The opposite side of the torque tube sits on two rollers. It was assumed, through experimental observations, that there was very little translational motion on the rollers along the torque tube axis. However, there was some rotation about the torque tube axis. This was captured in the model using a revolution joint (rotation only about torque tube axis). Fig.8 depicts the geometric locations of the three boundary conditions. The other components are connected to one another using a bonded contact condition when appropriate.

2.3.2 Wind load on the ANSYS model of heliostat

To validate the comparison between the simulated and the tested acceleration, one requirement should be met that the simulated wind load on the ANSYS heliostat model should be close to the real wind load. This paper tried to figure out the closest data of wind loads. In terms of wind load testing, the lowest anemometer was at 20 feet high, which was approximately the height of the top five facets of NSTTF heliostat. Therefore, the simulating data of wind velocity were picked out and compared with the actual wind velocity at 20 feet. Fig.9 shows the velocity comparison for facet 23 and facet 25. It can be seen that the simulated velocities agree with the testing velocity, which make the following comparison of the simulated acceleration and tested acceleration valid.

For the transient analysis, the time history of wind pressure was calculated according to the wind velocity time history and applied on the surface of 25 facets of the heliostat model. To reduce the computing time of transient analysis, the time length is two seconds. The wind load applied on the heliostat model begins with zero and goes up in slope curve during the first second as the pre-load. Then the simulated wind pressure is applied in the last second. The sampling frequency is 100Hz.

2.3.3 Modal-to-Transient analysis

Using the modal analysis, the mode shapes and the angular frequency of each mode were obtained such that the Rayleigh damp coefficient could be calculated and used for the transient analysis. Transient analysis was used to determine the dynamic response of the structure under the transient loads. The dynamic response includes deformation, stress, and strain and so on. In this paper, the directional deformation in the wind direction was obtained under the wind pressure time history. The transient dynamic analysis was based on the Equation (14). It can be seen from the equation that both inertial force and damping force are considered at any specified time. ANSYS programing solves the equation at discrete time points using time integration method.

[M]{u}+[C]{u}+[K]{u} = {F(t)} (14)

Fig.10 displays the Y directional deformation of the heliostat. The maximum deformation was about 6 mm which did not cause the facets to intersect. Fig. 10 is just a scaled deformation plot to show the deformation status visually. The Y axis is in the wind direction. The angle is 63.04 degrees between the surface of the heliostat and the wind direction. The elevation angle of the heliostat is 45 degrees. In addition, directional deformation time history of all points at which accelerometers were mounted was obtained and then converted to directional acceleration through twice differentiation using MATLAB.

3. Results and Discussion

Fig.10 Directional deformation of NSTTF heliostat

To compare the simulated acceleration with the tested acceleration, the acceleration time history was converted to the aptitude spectrum changing with the frequency. Fig. 11 displays the comparison of acceleration in frequency domain at some locations where accelerometers were mounted. Each sub-figure shows the dominant frequency which defines the frequency when the maximum amplitude of the acceleration occurs. It can be seen from the sub-figures that the shapes of both simulated and tested acceleration are similar, while the amplitude of the simulated acceleration is larger than that of the tested acceleration. The dominant frequencies at all mounting points for accelerometers were analysed and compared between the tested and simulated results. At the location of accelerometer No. 31, No.41 and No.202, the simulated and tested dominant frequencies are good matches. For the simulated acceleration, the dominant frequency value of 3.16 Hz had the highest probability of occurrence. For the tested acceleration, the corresponding frequency was 3.18 Hz.

The numerical simulation of wind loads on heliostats is the important part of work in this paper. To reduce the error of the simulation is of critical importance to get the correct dynamic response and then to obtain the good comparison with the tested results. During the simulating process, the wind load shape factor, the wind direction and the target spectrum may affect the simulation accuracy. The wind load shape factor is mainly related to the shape and dimension of the heliostat structure. Generally speaking, the accurate value of the factor can be obtained using wind tunnel test. In this paper, the factor was determined by experience value which may cause an error to some

И 0.04

Accelerometer--21Y Simulated Acceleration Tested Acceleration

Frequency (Hz)

a) Accelerometer No.21

0.15 -

3 0.05

Accelerometer--31Y Simulated Acceleration Tested Acceleration

10 20 30

Frequency (Hz)

b) Accelerometer No.31

Ф 0.2

Ф 0.1

s 0.02

Fig. 11 Comparison of acceleration in frequency domain

degree. In terms of the effect of wind direction on the simulation accuracy, the wind direction angle includes elevation angle and azimuth angle. It can be seen from the tested data that the elevation angle ranged from -20.05 to 33.68 degrees and the azimuth angle from 182.6 to 231.79 degrees. The average angles are 0.32 degree and 205.89 degrees respectively. In this paper, the elevation angle was ignored and only considered the average azimuth angle for the simulation of wind load. Besides, the choice of target spectrum is important for the simulation. As seen in Fig.4, Davenport spectrum is appropriate to be a target for the simulation. To obtain better simulation of wind load, the calculated wind profile around the heliostat should be accurate. Therefore, wind speed should be measured at more points, such as at the heights of 2m, 4m, 6m, 8m, 10m and 12m respectively, so that we can get the wind profile which is much closer to the actual situation.

Therefore, to analyze the dynamic response of the heliostat under wind load excitation, the wind load fluctuation analysis is firstly important. In future work, more wind load research methods should be explored to obtain the wind pressure time history such as CFD simulation, wind tunnel test or measurement on the real heliostat. For the method of numerical simulation described in this paper, the simulation accuracy should be further improved.

4. Conclusions

Numerical simulation methods for wind load and dynamic response analysis of heliostats have been investigated in this paper. The analysis result of the dynamic response was compared with the test result. The conclusions may be drawn as follows: 1) Auto-Regressive model can be used to simulate the wind load fluctuations on heliostats, which provides a method of wind load calculation for heliostat designers. 2) For the simulated acceleration which was obtained through transient analysis in ANSYS, the dominant frequency value of 3.16 Hz had the highest probability of occurrence. For the tested acceleration, the corresponding frequency was 3.18 Hz. This means the transient analysis using ANSYS is available to obtain dynamic response of heliostats. 3) Although the dominant frequency values of the simulated acceleration is close to that of the tested acceleration, it is still necessary to perform further research on improving the accuracy of simulating wind load fluctuations.

Acknowledgements

The authors would like to thank Sandia National Laboratories (SNL) for creating research conditions for this work and also thank Key Laboratory of Solar Thermal Energy and Photovoltaic System of Chinese Academy of Sciences for great supports. Special thanks will be given to the reviewers of this paper for their work and time.

References

[1] Bo Gong, Zhengnong Li, Zhifeng Wang, Yingge Wang. Wind-induced dynamic response of Heliostat. Renewable Energy 2012; 38:p206-213.

[2] Bo Gong, Zhifeng Wang, Zhengnong Li, Chuncheng Zang, Zhiyong Wu. Fluctuating wind pressure characteristics of heliostats. Renewable Energy 2013; 50:p307-316

[3] S. Huss, Y.D. Traeger, Z. Shvets, M. Rojansky, S. Stoyanoff, J. Garber, Evaluating effects of wind loads in heliostat design, SolarPACES, 2011

[4] C. K. Ho, D. T. Griffith, A. C. Moya, J. Sment, J. Christian, J. Yuan, and P. Hunter. Dynamic testing and analysis of heliostats to evaluate impacts of wind on optical performance and structural fatigue. Proceedings of SolarPACES 2012 International Conference, Marrakech, Morocco, September 11-14, 2012

[5] D. T. Griffith, A. C. Moya, C. K. Ho, P. S. Hunter, Structural dynamics testing and analysis for design evaluation and monitoring of heliostats, Proceedings of ASME 2011 5th International Conference on Energy Sustainability & 9th Fuel Cell Science, Engineering and Technology Conference, Washington, DC, USA , Aug.7-10, 2011

[6] D. Hiriart, J.L. Ochoa, and B. Garcia, Wind power spectrum measured at the San Pedro Martir Sierra, Revista Mexicana de Astronomiay Astrofisica,2001; 37, 213-220

[7] Linghua Zhang and Baoyu Zheng, Random signal processing, Tsinghua University Press, Beijing, 2003

[8] Nan Li, Study on numerical simulation of wind load and wind vibration control of large-span structures, Dissertation of graduate, Tianjin University, 2005