modelling of design of solar steam generating collector fields

Modelling and Design of Direct
Solar Steam Generating Collector
Fields
M. Eck
W.-D. Steinmann
German Aerospace Center (DLR),
Institute of Technical Thermodynamics,
Pfaffenwaldring 38-40, 70569
Stuttgart, Germany

The direct steam generation (DSG) is an attractive option regarding the economic improvement of parabolic trough technology for solar thermal electricity generation in the
¨
multi megawatt range. According to Price, H., Lupfert, E., Kearney, D., Zarza, E., Cohen,
G., Gee, R. Mahoney, R., 2002, ‘‘Advances in Parabolic Trough Solar Power Technology,’’
J. Sol. Energy Eng., 124 and Zarza, E., 2002, DISS Phase II-Final Project Report, EU
Project No. JOR3-CT 980277 a 10% reduction of the LEC is expected compared to conventional SEGS like parabolic trough power plants. The European DISS project has proven the feasibility of the DSG process under real solar conditions at pressures up to
100 bar and temperatures up to 400°C in more than 4000 operation hours (Eck, M.,
¨
Zarza, E., Eickhoff, M., Rheinlander, J., Valenzuela, L., 2003, ‘‘Applied Research Concerning the Direct Steam Generation in Parabolic Troughs,’’ Solar Energy 74, pp. 341–
351). In a next step the detailed engineering for a precommercial DSG solar thermal power plant will be performed. This detailed engineering of the collector field requires the consideration of the occurring thermohydraulic phenomena and their influence on the stability of the absorber tubes. ͓DOI: 10.1115/1.1849225͔

Introduction
The modelling of the DSG process in parabolic troughs is of special interest for the detailed design of the collector field. To identify critical process conditions that may cause e.g. an overheating of the absorber tubes a tool has been developed that considers all thermohydraulic aspects relevant for the design of a collector field such as the flow pattern in the evaporation section, the pressure loss and the heat transfer.
A design tool has been developed at DLR calculating all relevant process parameters including pressure drop, temperature field and stress in the absorber tubes. The models implemented in this design tool have been validated in detail at the DISS test facility under real DSG conditions for pressures between 30 and
100 bar and inner diameters between 50 and 85 mm. The models have been implemented into a MATLAB® program to allow for a first quick determination of critical process conditions. Once critical process conditions have been identified the FEM package
ANSYS® is used for a detailed investigation. This article summarizes the models used and shows the design procedure for a DSG collector field. The design program has proven to be a reliable tool for the detailed design of DSG collector fields.

heat transfer coefficient the knowledge of the flow pattern is essential for the calculation of the temperature distribution in the absorber tubes.
In practice various flow pattern maps are used for the determination of the present flow pattern. A frequently used flow pattern map is that of Taitel and Dukler ͓1͔. This flow pattern map is valid for horizontal and slightly inclined tubes. It has been developed for a two-phase flow without heat supply. To consider the effect of heating the models of Taitel and Dukler have been extended by some empirical equations ͓2͔. At DLR a different approach was chosen. During the experiments at the PRODISS and DISS test facility mainly wavy and annular flow occur ͑see Fig. 1͒. An empirical correlation has been developed that predicts the transition from wavy to annular flow. This model has been validated in a wide range of process parameters including different tube diameters ͓3͔. The mass flux density where the transition occurs can be calculated according to Eq. ͑1͒.

ͩ

˙ m g,t ϭ ͑ 46.6ϩ0.595pϩ0.0119p 2 ͒ 1ϩ1.3

˙ q 56

ͪ

(1)

Determination of the Boundary Conditions
A proper design of the collector field has to guarantee that the stress and the temperature in the absorber tubes do not exceed the limits of the absorber material used. To determine the stress and the temperature field the most important boundary conditions such as heat transfer coefficients, phase distribution, operation pressure, flux distribution and flux density have to be known. In the subsequent sections the models for the determination of these boundary conditions will be presented.
Flow Pattern. In the evaporator section of the collector loop a two-phase flow occurs. Depending on the process parameters different flow pattern and thus different phase distribution will occur. Since the present phase distribution determines the present
Contributed by the Solar Energy Division and presented at the ISEC2004 Portland, Oregon, July 11–14, 2004 of THE AMERICAN SOCIETY OF MECHANICAL
ENGINEERS. Manuscript received by the ASME Solar Division; April 27, 2004; final revision August 10, 2004. Associate Editor: J. Davidson.

Journal of Solar Energy Engineering

Equation ͑1͒ is an equation between quantities where the pressure p is given in bar, the heat flux density q is given in kW/m2.
The mass flux density m t is calculated in kg/m2 s.
Wetting Angle. In case of annular flow the complete inner circumference of the absorber tube is wetted by saturated water whereas at wavy flow the inner circumference is only partly wetted. In this case the boundary between the wetted and the unwetted region is of special interested. This boundary is expressed by the wetting angle ␸ wet ͑see Fig. 2͒. To determine the wetting angle of a wavy flow the liquid level for a fictitious stratified flow is calculated first. Afterwards this liquid level is multiplied by a so called wave factor to get the wetting level of a wavy flow.
At stratified flow the liquid level depends on the void fraction ␧
͑the ratio of the cross section occupied by steam to the total inner cross section of the tube͒ that is determined by the equation of
Rouhani ͓4͔.

Copyright © 2005 by ASME

AUGUST 2005, Vol. 127 Õ 371

Downloaded From: http://solarenergyengineering.asmedigitalcollection.asme.org/ on 07/14/2013 Terms of Use: http://asme.org/terms

Fig. 1 Main flow patterns for direct steam generation

␧ϭ

ͫ

ͩ

˙
˙
˙ x x
1Ϫx
˙ ϩ ͑ 1ϩ0.12͑ 1Ϫx ͒͒
␳g
␳g
␳l
ϩ

˙
1.18͑ 1Ϫx ͒͑ g ␴ ͑ ␳ l Ϫ ␳ g ͒͒ 0.25
˙
m ␳ l0.5

ͬ

ͪ

Ϫ1

(2)

Fig. 3 March of the specific pressure loss along the collector loop „d i Ä50 mm, MÄ1 kgÕs…

The normalized liquid level h l /d i is determined by the implicit equation ͑3͒.

ͫ

ͩ

ͩ

hl sin 2 arccos 2 Ϫ1 di 1 hl ␧ϭ arccos 2 Ϫ1 Ϫ
␲
di
2

ͩ

ͪ

ͪͪ ͬ

(3)

The wetting angle finally results from Eq. ͑4͒.

ͩ

hl
Ϫ1
␸ wetϭarccos d i /2

ͪ

ͩ ͪ ͩ ͪ dp dl

h wϭ f h l

˙
˙
RϭAϩ3.43x 0.685͑ 1Ϫx ͒ 0.24

ͩ

ϫ 1Ϫ

ͩ ͪ

1 ␳ 2 w , ϭ␭ d 2
1ph

␩g
␩l

ͪ

0.89

Ϫ0.25

␭ϭ0.316 Re

with
˙
˙
Aϭ ͑ 1Ϫx ͒ 2 ϩx 2
Frl ϭ

(6)

According to the investigations during the European DISS project ͓5͔ the pressure loss of the occurring two phase flow is determined best with the correlations of Friedel. For a two-phase flow the pressure drop is defined as the product of the singlephase water flow and the two-phase flow multiplier R.

(7)
1ph,l

ͩ ͪ ͩ ͪ
␯g
␯l

0.8

•

␩g
␩l

0.22

FrlϪ0.47 Wel0.0334

(5)

Pressure Drop. In the turbulent single phase water or steam flow the pressure drop is determined by the Blasius equation ͑6͒. dp dl

2ph

dp dl The single phase water flow is calculated according to Eq. ͑6͒ assuming that the complete mass flux is liquid. The two-phase flow multiplier is calculated with Eq. ͑8͒:

(4)

At a wavy flow the wetting angle depends on the wave height h w . The ratio of wave height to the liquid level h l is expressed by the wave factor f which is in the range of 1.7 to 2.

ϭR

ͩ ͪ
␯ g␰ g
␯ l␰ l

(8)

0.8

(9)

˙
͑ m␯l͒2 gd i

(10)

˙ m 2d i␯ l
␴

(11)

Wel ϭ

The parameters ␰ g and ␰ l are calculated for turbulent steam flow (RegϾ1055) as follows:

ͫ

␰ g ϭ 0.86859 log

ͩ

Reg
1.964 log͑ Reg ͒ Ϫ3.8215

ͪͬ

Ϫ2

(12)

For laminar steam flow (RegϽ1055)

␰ gϭ

64
Reg

For turbulent water flow (RelϾ1055)

ͫ

␰ l ϭ 0.86859 log

ͩ

Rel
1.964 log͑ Rel ͒ Ϫ3.8215

(13)

ͪͬ

Ϫ2

(14)

For laminar water flow (RelϽ1055)

␰ lϭ

Fig. 2 Schematic cross section of an evaporation tube

372 Õ Vol. 127, AUGUST 2005

64
Rel

(15)

Figure 3 displays the local specific pressure loss at different axial positions of a collector loop for different operation pressures. Transactions of the ASME

Downloaded From: http://solarenergyengineering.asmedigitalcollection.asme.org/ on 07/14/2013 Terms of Use: http://asme.org/terms

Fig. 4 Heat transfer coefficient as a function of the spec. enthalpy for different heat flux densities. In two-phase region only the heat transfer coefficient in the wetted and heated region is
˙
displayed. „pÄ100 bar, M Ä1 kgÕs, d i Ä50 mm…

Fig. 5 Heat transfer coefficient as a function of the spec. Enthalpy for different pressures. In the two-phase region only the heat transfer coefficient in the wetted and heated region is dis˙
˙
played. „q Ä40 kWÕm2 , M Ä1 kgÕs, dÄ50 mm…

The two-phase flow pressure loss in the bends is determined by the equation of Chisholm ͓6͔ calculating the friction factors according to ͓7͔.

cient as a function of the spec. Enthalpy along the collector loop for different boundary conditions.
At this point all boundary conditions at the inner circumference are defined as function of the known process parameters, the operation pressure, and the fluid mass flux. In the following sections the boundary condition at the outer circumference has to be defined. Heat Transfer Coefficient. The heat transfer coefficient for the single-phase water or steam flow is calculated according to the equation of Dittus–Boelter ͓8͔.
␭
␣ 1phϭ0.0235 Re0.8 Pr0.48 d (16)

In case of the two-phase water-steam flow three different heat transfer phenomena with different heat transfer coefficients can be identified. 1. Heat transfer to the steam in the unwetted region;
2. Heat transfer to the water in the wetted and heated region;
3. Heat transfer to the water in the wetted and unheated region.

Collector Efficiency. Parabolic trough collectors only use the direct normal irradiation ͑DNI͒. The DNI is reflected on the absorber tube. Due to optical losses at the reflector, the glass envelope and the absorber surface only a certain fraction of the DNI is absorbed by the absorber tube. In addition the absorber tube has thermal losses due to convection and radiation. All these losses are determined by the collector efficiency. Once the collector efficiency is known the heat flux on the outer surface of the absorber tube is known.

The heat transfer coefficient for the steam phase and the water phase in the unheated region can be determined with Eq. ͑16͒. In the wetted and heated region nucleate boiling occurs. According to Gungor and Winterton ͓9͔ the heat transfer coefficient ␣ 2ph,l is determined by a first term taking into account the convective heat transfer ␣ conv and a second term for the heat transfer by the nucleate boiling ␣ NB .

␣ 2ph,l ϭE ␣ convϩS ␣ NB

(17)

The supression factor S is determined as follows:
Sϭ ͓ 1ϩ1.15•10Ϫ6 E 2 Ret1.17͔ Ϫ1

(18)

The enhancement factor E is calculated according to Eq. ͑19͒.
Eϭ1ϩ24000Bo 1.16ϩ1.37X Ϫ0.86 a (19)

With the boiling number Bo.
Boϭ

˙ q ˙ m tot͑ h Љ Ϫh Ј ͒

(20)

˙ where m tot is the total mass flux density. The Martinelli parameter
˙
X a is a function of the fluid properties and the steam quality x .
The convective heat transfer coefficient again is calculated according to Eq. ͑16͒. The heat transfer coefficient due to nucleate boiling is determined according to Cooper ͓10͔.
Figures 4 and 5 display the course of the heat transfer coeffiJournal of Solar Energy Engineering

Fig. 6 Collector efficiency as a function of the difference between operation and ambient temperature „TÀT a … for different values of the DNI „␸Ä0 deg, Cermet with Vacuum, LS-2 collector…

AUGUST 2005, Vol. 127 Õ 373

Downloaded From: http://solarenergyengineering.asmedigitalcollection.asme.org/ on 07/14/2013 Terms of Use: http://asme.org/terms

not be higher than 50 K above the fluid temperature ͓16͔. To allow for a very fast and reliable determination of the highest absorber temperature a simplified analytical solution for the calculation of the circumferential temperature distribution on the medium center line has been developed and implemented into the design tool.
Based on this result the maximum and minimum temperature at the outer surface will be calculated.
Temperature on the Medium Center Line. The analytical calculation of the temperature is based on a division of the cross section in four segments ͑see Fig. 8͒.
1.
2.
3.
4.

Fig. 7 Typical heat flux distribution along the outer surface of an absorber tube and its approximation by a gaussian and a rectangular distribution with ␴Ä60 deg

Empirical correlations describing the collector efficiency for several collector models exist. For LS-2 collectors the efficiency has been investigated in ͓11͔ where different equations have been derived for different absorber configurations. Formulations for the
EuroTrough collector are presented in ͓12͔.
The typical course of the efficiency of a parabolic trough collector is displayed in Fig. 6 for an LS-2 collector. The efficiency of the collector decreases with an increasing operation temperature and a decreasing DNI. The maximum efficiency at TϪT a ϭ0, where no thermal losses occur, is called optical efficiency.
These functions are valid for an incidence angle of 0 deg ͑the sun is perpendicular to the aperture area of the collector͒. The effect of different incidence angles is considered by the incidence angle modifier ͑IAM͒. Again empirical correlations have been developed for different collectors. The IAM for an LS-2 collector is given in ͓11͔. Functions for an LS-3 collector are given in ͓13͔ and ͓14͔.
Every parabolic trough collector mentioned above can be used for direct steam generation. For the determination of the heat input into the absorber tube the designer of the collector field has to choose the efficiency and IAM curve for the collector used.

Wetted and heated;
Unwetted and heated;
Unwetted and unheated;
Wetted and unheated.

In case of the preheating and the superheating section with the according single-phase flow only two different sections can be identified, which are particular cases of the common case.
It is assumed that the boundary conditions in the different sections are known and constant. For the exact segmentation into the four sections the wetting angle and the position of the rectangular heat flux distribution must be known. The position of the heat flux distribution depends on the position of the sun that can be determined according to e.g. ͓17͔.
Figure 9 is used for the derivation of the accounting equations.
At the outer surface ⌬Q sol represents the heat input by the solar energy, ⌬Q conv represents the heat transfer at the inner surface and
Q ␭ resp. ⌬Q ␭ the conductive heat transfer to adjacent segments. It can be written that
˙
˙
˙
dQ ␭ ϩdQ convϭdQ sol

(21)

with
2

d T
˙
dQ ␭ ϭϪ␭s⌬z 2 dy dy (22)

˙ dQ convϭa⌬z ͑ TϪT f ͒ dy

(23)

˙
˙
dQ solϭ⌬zq soldy

(24)

In the heated region the heat flux density is:
˙
q sol,h ϭ ␩ col DNI d ap

360 d m␲͑ ␸ eϪ ␸ b ͒

(25)

Flux Distribution on the Outer Surface. Once the heat flux on the outer surface is known its distribution has to be determined.
A typical flux distribution at the outer circumference of an absorber tube for parabolic trough collectors is displayed in Fig. 7
͑taken from ͓15͔͒. This distribution shows a local minimum at an angle of 180 deg caused by the shadow of the absorber tube on the mirror facets. This distribution is similar to a Gaussian distribution with a standard deviation of 60 deg as displayed in Fig. 7.
Additionally the approximation by a rectangular distribution with its edges at 100 deg and 260 deg is presented. In all cases the integral of the curve is the same.
FEM analyses of the resulting temperature distribution in the absorber cross section have shown that there is no significant difference between the temperature distribution calculated with a flux distribution according to the Gaussian distribution and the one calculated with the typical parabolic trough distribution. Accordingly the Gaussian distribution is chosen for the accurate
FEM analysis whereas the rectangular distribution is used for a simplified analytical solution presented in the next section.

Temperature Distribution
The maximum temperature of the absorber tube is the most critical parameter of a DSG collector loop. This temperature may
374 Õ Vol. 127, AUGUST 2005

Fig. 8 Schematic illustration of an absorber cross section with the four different sections

Transactions of the ASME

Downloaded From: http://solarenergyengineering.asmedigitalcollection.asme.org/ on 07/14/2013 Terms of Use: http://asme.org/terms

Fig. 9 Developed view of an absorber tube segment for the derivation of the analytical solution. The depth of the Segment is ⌬z.

It is assumed that the thermal radiation at the outer surface is constant and equivalent to:
˙
q sol,unhϭϪ

͑ ␩ optϪ ␩ ͒ DNId ap d m␲

(26)

The differential equation for the determination of the temperature distribution is given by: d 2T dy 2

ϪK 2 TϭϪ

1
˙
q ϪK 2 T f
␭s sol

(27)

with
K 2ϭ

␣
␭s

(28)

And the solution:
T ͑ y ͒ ϭC 1 cosh͑ Ky ͒ ϩC 2 sinh͑ Ky ͒ ϩ

˙ q sol ϩT f
␣

(29)

To determine the two coefficients C 1 and C 2 two boundary conditions per segment have to be known. Since the four sections are connected to each other, the boundary conditions at the end of one section is the same as at the beginning of the next section. At this location the according temperatures as well as their gradient have to be the same. For example:
T 1 ͑ yϭl 1 ͒ ϭT 2 ͑ yϭ0 ͒

(30)

dT 1 ͑ yϭl 1 ͒ T 2 ͑ yϭ0 ͒ ϭ dy dy (31)

and

a temperature of 310°C whereas the steam is slightly superheated to a temperature of 320°C. The circumferential angle ␸ is counted clockwise from the crest.
The marches of the temperature are very similar in Fig. 10. The maximum deviation between the two different calculation procedures is approx. 5 K. For a first rough estimation of the temperature field this deviation is acceptable.
Maximum Temperature. To assess the thermal load of the absorber tube not the temperature distribution on the medium center line is of importance but the maximum temperature at the outer surface of the absorber tube, which is the hottest point, that can be determined from the known temperature on the medium center line. Therefore it is assumed that at the location with the maximum temperature the temperature gradient in circumferential direction is zero. Thus the heat flow in circumferential direction is zero. Accordingly only a heat flux in radial direction occurs. The maximum temperature is determined with the equation for the steady-state heat conduction in radial direction. Accordingly the maximum temperature at the outer surface is calculated by: ro rm
T max,o ϭT max,m ϩ
͑ T max,m ϪT max,i ͒ rm ln ri ln

The gradient of the temperature is defined by: dT 1 ͑ y ͒ ϭC 1 K sinh͑ Ky ͒ ϩC 2 K cosh͑ Ky ͒ dy Fig. 10 Comparison of the temperature along the medium center line of the absorber cross section calculated with the FEM package ANSYS® and the analytical solution. „two-phase-flow, heated from the side…

(32)

Accordingly there are eight equations for the determination of the eight coefficients C i . It is possible to simplify the equations by considering the special characteristics of the hyperbola function: sinh͑ 0 ͒ ϭ0

(34)

With the help of this simplified procedure a first rough assessment of a DSG collector loop is possible and critical process situations can be identified. In a subsequent step critical situations will be investigated in more detail using the FEM package
ANSYS® .

(33)

cosh͑ 0 ͒ ϭ1

(35)

The values for the eight coefficients C i are found by solving the set of eight linear equations. Knowing the coefficients C i for the different segments, it is possible to calculate the temperature on the medium center line. By connecting the four segments the march of the temperature along the medium center line of the investigated cross section is determined.
Figure 10 displays the march of the absorber temperature calculated with the FEM package ANSYS® and the analytical solution. For the calculation it was assumed that the liquid phase has
Journal of Solar Energy Engineering

Design Tool
The models presented have been implemented in a simulation tool using the programming environment MATLAB® allowing for a fast investigation of different collector loop configurations. Before starting the calculation run the boundary conditions such as pressure and temperature at the loop outlet, the direct normal irradiation, the position of the sun and the geometry of the collector loop have to be defined. In case of the recirculation mode the recirculation rate and the use of an injection cooler in the superheating section has to be defined too. The once-through mode
AUGUST 2005, Vol. 127 Õ 375

Downloaded From: http://solarenergyengineering.asmedigitalcollection.asme.org/ on 07/14/2013 Terms of Use: http://asme.org/terms

Table 1

Parameters of the ET-II collector

Parabola Width ͓m͔
Overall Length of a single Collector ͓m͔
Outer Diameter ͓mm͔
Inner Diameter ͓mm͔
Length of connecting pipe between two adjacent collectors ͓m͔
Number of 90° elbows in U-bends ͓°C͔
Peak Optical Efficiency ͓%͔
Roughness Factor of the Absorber Pipes ͓m͔

5,76
98,5
70
55
4
8
76,5
4EϪ05

Table 2 Main Parameters of the Solar Field
Length of Evaporator ͓m͔
Length of Superheater ͓m͔
Outlet Pressure ͓bar͔
Outlet Temperature ͓°C͔
Electric Power ͓MW͔
Number of parallel rows

Fig. 11 Pressure along the collector loop for the different operation modes „pÄ100 bar, T out Ä400°C, d i Ä50 mm, DNIÄ800
WÕm2…

requires the definition of the inlet temperature and the injection mode the definition of the number and distribution of injection cooler along the collector loop.
For the calculation the collector loop is divided into a number of segments with a step size of usually 5 m. Starting from the end of the collector loop, the collector efficiency, the specific enthalpy, the pressure loss and the fluid temperature for all investigated collector positions are calculated. Based on the calculated fluid properties and flow conditions the distribution of the heat transfer coefficient and the temperature of the absorber tube along the collector loop is determined. In addition it is possible to display the evaporation path in the Taitel–Dukler diagram.
As an example Fig. 11 displays the pressure and Fig. 12 the maximum temperature difference along the collector loop for the three different operation modes. According to Fig. 11 the total pressure drop is less than 2 bar for a length of the collector loop of
500 m and an outlet pressure of 100 bar.
As shown in Fig. 12 in all cases the highest temperatures do not occur in the evaporation but in the superheating section. This is

800
200
70
410
5
7

caused by the occurrence of an annular flow in the evaporation section with the according good cooling of the absorber wall.
Lay-Out of a DSG Collector Loop. In the previous sections the fundamentals of DLR’s design tool are presented. This tool is used during the design phase of DSG collector loops. ͑Since parabolic trough collectors operated with a heat transfer fluid such as oil are similar to the preheating section of a DSG collector loop,
DLR’s design tool is also applicable to oil collector loops.͒ The detailed design of collector fields also requires the investigation of the absorber tube and of the interconnecting piping and the necessary fittings in more detail. For the detailed investigation of the absorber tube the FEM package ANSYS® and for the piping and fittings and the overall collector field performance the heat balance calculation program IPSEpro® is used. The use of these three programs for the investigation of DSG collector fields is presented in the following section.
The investigated collector field is operated in recirculation mode as proposed in ͓18͔. To guarantee a sufficient cooling of the superheating collectors a length of 200 m is chosen fo