Numerical investigation of turbulent forced convection heat transfer in pillow plates

Numerical investigation of turbulent forced convection heat transfer in pillow plates

International Journal of Heat and Mass Transfer 94 (2016) 516–527 Contents lists available at ScienceDirect International Journal of Heat and Mass T...

3MB Sizes 664 Downloads 895 Views

International Journal of Heat and Mass Transfer 94 (2016) 516–527

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Numerical investigation of turbulent forced convection heat transfer in pillow plates M. Piper a, A. Zibart a, J.M. Tran a, E.Y. Kenig a,b,⇑ a b

Chair of Fluid Process Engineering, Faculty of Mechanical Engineering, University of Paderborn, Paderborn, Germany Gubkin Russian State University of Oil and Gas, Moscow, Russian Federation

a r t i c l e

i n f o

Article history: Received 3 August 2015 Received in revised form 27 October 2015 Accepted 5 November 2015 Available online 17 December 2015 Keywords: Pillow-plate Fluid dynamics Heat transfer CFD

a b s t r a c t Pillow-plate heat exchangers (PPHE) are characterized by the wavy ‘‘pillow-like” surface. Being a promising alternative to conventional equipment, they have recently received increased attention in the process industry. However, the knowledge on the flow and heat transfer characteristics in pillow-plate channels is scarce. Such knowledge is necessary for the development of successful design principles for PPHE. In this study, turbulent single-phase flow in pillow plates was investigated using CFD. The results for pressure loss and heat transfer coefficients were successfully validated against our own experiments. The simulations show that the flow in pillow plates is characterized by two distinct regions, namely a meandering core, which dominates heat transfer, and recirculation zones, which occur in the wake of the welding spots. The influence of the pillow-plate geometry on their thermo-hydraulic performance was also investigated, and the most favorable designs were identified. Based on the simulated flow field simple geometrical improvements are proposed, which significantly improve the thermo-hydraulic performance of pillow plates. Ó 2015 Elsevier Ltd. All rights reserved.

1. Introduction The development of successful design methods for pillow-plate heat exchangers (PPHE) requires knowledge of the flow and heat transfer characteristics in the pillow-plate channels. Whereas such knowledge is largely available for conventional equipment (e.g. shell-and-tube and corrugated-plate heat exchangers), it is very scarce for PPHE. The lack of publicly available, reliable design methods for PPHE and the absence of reference applications in industry hinders their widespread application. Fig. 1 shows a representative element of a pillow plate, characterized by its wavy ‘‘pillow-like” surface. The geometric parameters shown in Fig. 1 are characteristic for pillow plates. These are: the longitudinal (with respect to the flow direction) welding spot pitch sL , the transversal welding spot pitch sT , the diameter of the welding spots dSP , the maximal inner inflation height di and the thickness of metal sheets dP . The manufacturing of pillow plates involves the following steps: first two metal sheets are spot welded to each other over their entire surface using resistance or laser welding. The regular welding spot pattern is chosen ⇑ Corresponding author at: Chair of Fluid Process Engineering, Faculty of Mechanical Engineering, University of Paderborn, Paderborn, Germany. Tel.: +49 (0) 5251602408; fax: +49 (0) 5251602183. E-mail address: [email protected] (E.Y. Kenig). http://dx.doi.org/10.1016/j.ijheatmasstransfer.2015.11.014 0017-9310/Ó 2015 Elsevier Ltd. All rights reserved.

according to requirements on thermo-hydraulic performance and structural stability. Second, the edges of the pillow plates except for the connecting ports are welded by roll seam welding. Finally, the metal sheets are inflated by hydroforming, thus giving the pillow plate its characteristic pillow-like surface. The welding spots and the channel waviness promote lateral mixing, which enhances heat transfer and results in good thermo-hydraulic performance of pillow plates. Because of their fully welded construction, pillow plates are hermetically sealed and offer high structural stability. The light weight and compactness of PPHE make them ideal as top condensers in distillation columns, where they can be implemented directly into the column head, thus reducing material costs (no need for extra housing) and required space. Literature on pillow plates is scarce. Few published works dealing with single-phase flow in pillow plates are by Mitrovic and Peterson [1], Mitrovic and Maletic [2], Piper et al. [3] and Tran et al. [4,5]. Studies [1,4,5] are of an experimental and [2,3] of a numerical nature. The latter papers focused on gaining a deeper understanding of fluid dynamics and heat transfer in pillow-plate channels. In [2], a comprehensive CFD-study of the flow in pillow plates over the range 50 6 Re 6 3800 was accomplished. A laminar model was used in all simulations, although the flow in the middle and upper Reynolds number range was turbulent. As denoted by Maletic [6], this simplification resulted in an underestimation of

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

517

Nomenclature A bD d h k lchar lR Nu p Pr Q_ q_ Re s sr T u us V_ _ W x; y; z ~ y yþ

area free diagonal distance between welding spots diameter heat transfer coefficient turbulent kinetic energy characteristic length free longitudinal distance between welding spots Nusselt number pressure Prandtl number heat flowrate heat flux Reynolds number welding spot pitch reduced welding spot pitch temperature velocity shear stress velocity volumetric flowrate pumping power Cartesian coordinates coordinate normal to the wall dimensionless coordinate normal to the wall

Greek symbols Dp pressure loss DT temperature difference d channel height e thermo-hydraulic efficiency e normalized thermo-hydraulic efficiency  turbulent dissipation rate f friction coefficient k thermal conductivity m kinematic viscosity q density

pressure loss and especially heat transfer rate, compared to the experimental results of [1]. Maletic [6] also performed simulations with the standard k— turbulence model; however, this led to largely overestimated heat transfer coefficients. This shows that the use of an appropriate turbulence model is necessary. A further simplification employed in [2] was the approximation of the wavy pillow-plate surface by some trigonometric function. This resulted in cross-sectional areas differing from reality (see [7]).

sw w

x

wall shear stress dimensionless ‘‘two-zone” parameter specific dissipation

Subscripts A area B bulk C circular cs cross-section E elliptic char characteristic D drag exp experiment d thickness of viscous sub-layer h hydraulic i inner k index L longitudinal m mean mc meandering core P plate (metal sheet) Q heat R friction Re Reynolds ref reference RZ recirculation zone SP welding spot sim simulation sym symmetry T transverse # temperature tot total u velocity w wall Dp pressure loss

In [3], turbulent single-phase flow and heat transfer in pillow plates were studied over the range 1000 6 Re 6 8000 using CFD. The complex pillow-plate geometry was accurately reconstructed using forming simulations, which imitated the manufacturing process of pillow plates. The simulation results for pressure loss were validated against measurements of our group and showed a mean deviation of less than 5%. A validation of heat transfer coefficients gained by the CFD model used in [3] was still necessary. Furthermore, the influence of different pillow-plate geometries on the thermo-hydraulic performance was not investigated in [3]. In this work, a comprehensive CFD study on turbulent singlephase flow and heat transfer in pillow plates was carried out. The results gained in [3] were extended to a large variety of different pillow-plate geometries, thus gaining a more general understanding of the flow and heat transfer characteristics in pillow plates. Furthermore, the influence of the pillow-plate geometry on the thermo-hydraulic performance was investigated revealing the most favorable geometries. In addition, methods for optimizing the pillow-plate geometry with respect to improved thermohydraulic performance were proposed.

2. Geometry generation and summary of studied pillow plates

Fig. 1. Characteristic geometrical parameters of pillow plates. The realistic geometry was generated using forming simulations.

For the CFD simulations, a digital image of the pillow-plate geometry that defines the boundaries of the simulation domain is required. The prerequisite for the realistic description of the fluid

518

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

dynamics in pillow plates is an accurate reconstruction of the pillow-plate channel. In this work, the pillow-plate geometry was generated by forming simulations. This method, which is described in detail in [7], allows a realistic reconstruction of the pillow-plate geometry, with a deviation from the real pillow-plate surface of less than 4%. The pillow plates generated and investigated in this study are summarized in Table 1. The welding spot pattern can be characterized in dimensionless form by the following parameter:

sr ¼

2sL  dSP sT  dSP

ð1Þ

This parameter is equal to one for an equidistant welding spot pattern (2sL  dSP ¼ sT  dSP ), smaller than one for a ‘‘transversal” pattern (2sL  dSP < sT  dSP ) and greater than one for a ‘‘longitudinal” pattern (2sL  dSP > sT  dSP ). The hydraulic diameter dh is determined according to the method proposed in [7]. The pillow-plate 72/42/3/10, used in our experimental set-up (see [5]), is our reference geometry in this study. It has a typical triangular welding spot arrangement with arctanð2sL =sT Þ  60 , which is a pattern commonly found in industry. All other values of sT , sL and di presented in Table 1 are systematic variations of 72/42/3/10, e.g. sT and sL are interchanged or the inflation height di is doubled. 3. CFD-simulation The CFD simulations were performed using the commercial solver STAR-CCM+ by CD-Adapco, which is based on the finite

Table 1 List of investigated pillow-plate geometries. 2sL (mm)

sT (mm)

di (mm)

dSP (mm)

dh (mm)

sr ()

Abbreviation

42 42 42 42 42 42 50 57 60 60 64 64 68 68 72

72 72 57 57 42 42 42 42 42 42 42 42 42 42 42

3 6 3 6 3 6 6 6 3 6 3 6 3 6 3

10 10 10 10 10 10 10 10 10 10 10 10 10 10 10

4.43 8.57 3.84 7.31 3.37 6.32 6.72 7.31 3.99 7.52 4.24 8.07 4.38 8.36 4.43

0.52 0.52 0.68 0.68 1.0 1.0 1.25 1.47 1.56 1.56 1.69 1.69 1.81 1.81 1.94

42/72/3/10 42/72/6/10 42/57/3/10 42/57/6/10 42/42/3/10 42/42/6/10 50/42/6/10 57/42/6/10 60/42/3/10 60/42/6/10 64/42/3/10 64/42/6/10 68/42/3/10 68/42/6/10 72/42/3/10

volume method (FVM). The computational domain illustrated in Fig. 2 was obtained with the aid of forming simulations. The computational effort was reduced by exploiting all existing symmetries of the flow (see Fig. 2). A similar domain was used in [2]; however, we further reduced the domain by utilizing the flow periodicity. This was accomplished by applying periodic boundary conditions at the inlet and outlet boundaries of the channel. In this way, the velocity field is repeated after twice the longitudinal pitch (2sL ), which offers the advantage of limiting the investigation to the (periodically) hydrodynamically developed region. In industrialscale pillow-plate heat exchangers, such regions usually occupy a major part of the channels. Boundary effects and the influence of inlet and outlet including eventual distribution areas are aspects to be considered in future work. Thermally developed flow is characterized by a constant heat transfer coefficient in flow direction (h ¼ f ðyÞ ¼ const). In STARCCM+, this is achieved by scaling the self-similar temperature profiles at the inlet and outlet boundaries (only), in such a way that the heat transfer coefficient at these boundaries is the same. At the walls, we applied a no-slip boundary condition (u ¼ 0) and assumed constant wall temperature (T w ¼ const). The latter appears reasonable, because a major application area of PPHE is condensation, e.g. in top condensers in distillation columns, where this boundary condition is approximately fulfilled. The flow considered was single-phase, incompressible, steadystate, three-dimensional and turbulent, with constant physical properties. Turbulence was described statistically by the Reynolds-Averaged-Navier–Stokes (RANS) equations. The Reynolds stresses were computed using the realizable k— model (Shih et al. [8]) available in STAR-CCM+. In contrast to the standard k—, the realizable k— model is free of such problems as round-jet and stagnation point anomalies [9]. Furthermore, the realizable k— model shows superior performance in flows involving rotation, boundary layers under strong adverse pressure gradients, separation and recirculation [8]. It has been applied successfully to complex geometries, such as structured packings (e.g. [10,11]) as well as plate [12] and shell-and-tube heat exchangers [13]. Pressure loss and heat transfer coefficients obtained using this model agreed well with our experiments. Since flow separation was expected to occur in the pillow-plate channel, it was necessary to resolve the boundary layers appropriately. This was accomplished by using a two-layer formulation (see [14,15]), which allows the application of the turbulence model also in the viscosity dominated region. The boundary layer was decomposed into two regions: a turbulent and a viscous one. In the turbulent region, the realizable k— model is retained. In the viscous region, however, the turbulent dissipation rate and turbulent viscosity are functions of the distance normal to the wall. In this near-wall region, the one-equation model of Wolfshtein [16]

Fig. 2. Periodic computational domain of the pillow-plate channel with staggered welding spots.

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

is employed. Hence, only one transport equation is solved, namely the equation for the turbulent kinetic energy, while the dissipation rate  is evaluated by an algebraic expression. Finally, the values of the dissipation rate determined in the near-wall region are smoothly coupled with those computed far from the wall. In [3] we used a block structured grid for the simulations. The boundary layers were resolved using prismatic cells, which aligned tangentially with the wavy surface. In the bulk flow, a Cartesian grid was applied. Gradients in the directions x and y, i.e. other than normal to the wall, were resolved using a flow-oriented local grid refinement. With the release of a newer version of STAR-CCM+ (10.02.010), with improved meshing capabilities, we were able to move to structured body-fitted grids. This grid type offers some advantages over the grid used in [3], namely, more convenient and faster grid generation, reduced computational effort due to a lower number of cells needed, higher cell quality and improved convergence. It is worth noting that the simulation results obtained by these two grids were identical. An example of the body-fitted grid used in this work is illustrated in Fig. 3, showing how the cells align with the geometry contours, e.g. in the vicinity of the welding spots in the x—y-plane and close to the wavy channel wall. In addition to the testing of different grid types (as mentioned above), we also performed a grid-independence study in order to ensure an adequate grid resolution. This study was performed always for the highest Reynolds number (Re ¼ 8000) in each geometry. It was assumed that the grid resolution sufficient for Re ¼ 8000 would also be sufficient for the lower Reynolds numbers. Depending on the geometry, the total cell number could reach 15  106 cells. For properly resolving the boundary layers, the following steps were undertaken when constructing the grid:  The distance between the wall and the grid node closest to it was estimated using the dimensionless wall coordinate of the inner region of the turbulent boundary layer:

yþ ¼

~ us y

m

rffiffiffiffiffiffi with us ¼

sw q

ð2Þ

~ the coordinate norIn Eq. (2), us denotes the shear velocity and y mal to the wall. The grid node adjacent to the wall had a value yþ < 1.

519

 At least 3–4 cells are present in the viscous sub-layer (0 6 yþ 6 5), in order to resolve this region properly, ensuring a smooth transition to the buffer layer. ~k ) increased linearly (‘‘stretched”) in direc The cell thickness (Dy ~kþ1 =Dy ~k ) was tion normal to the wall. The stretching factor (¼ Dy set to a value of about 1.1. In the simulations of heat transfer with Pr  1, the conduction layer is thinner than the viscous sub-layer. Hence, a smaller yþ value for the grid node adjacent to the wall was used. For large Prandtl numbers, the thickness of the conduction layer (yþ #;d ) was approximately determined by the following relation [17]:

yþ#;d ¼ 15Pr1=3

ð3Þ

4. Process parameters definitions For the evaluation of the thermo-hydraulic characteristics of the flow in pillow plates, several process parameters are introduced. The mean Reynolds number of the flow in the pillow-plate channel is defined using the following expression:

uchar lchar

Re ¼

m



um dh

m

ð4Þ

The characteristic velocity uchar is represented by the xcomponent of the mean velocity in the channel. As the characteristic length lchar , the hydraulic diameter of the channel determined by the methods proposed in [7] is taken. The pressure loss coefficient is evaluated by the Darcy–Weisbach equation:

fDp ¼

2dh Dp

ð5Þ

qu2m 2sL

Pressure loss is calculated from the difference between the surface-averaged pressure values at the inlet and outlet boundaries. The hydraulic Fanning factor is determined using the surface averaged wall shear stress:

fR ¼

8sw qu2m

ð6Þ

The form drag coefficient is obtained by subtracting Eq. (6) from Eq. (5):

fD ¼ fDp  fR

ð7Þ

The Nusselt number is defined by

Nu ¼

hm dh k

ð8Þ

where hm represents the surface-averaged heat transfer coefficient:

hm ¼

q_ m Tw  Tm

ð9Þ

In Eq. (9), q_ m is the heat flux in the direction normal to the wall and T m is the adiabatic mean temperature in the channel evaluated under the assumption of constant physical properties as follows:

Z

uT dA T m ¼ ZAcs

:

ð10Þ

u dA Acs

Fig. 3. Example of structured body-fitted grid used in the simulations. The grid here is extra coarse for illustration purposes.

Here Acs denotes the cross-sectional area. In order to consistently compare computed and measured heat transfer coefficients, we used T m instead of the temperature T sym at the channel axis (symmetry plane ðx; y; 0Þ), since usually only T m is experimentally accessible, especially in complex geometries.

520

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

between the two ports was equal to 8 times the longitudinal welding spot pitch (i.e. 288 mm). The pressure loss measurement was performed under isothermal conditions using water at T ¼ 25  C as a working fluid. The pillow-plate used in the experiments has industrially relevant characteristic dimensions: sT ¼ 42 mm, 2sL ¼ 72 mm, dSP ¼ 10 mm and di ¼ 3 mm. It was manufactured from austenitic stainless steel plates of material 1.4541 (AISI 321), with a surface finishing of quality 2B (DIN EN 10088-2), which has a typical mean roughness Ra  0:1—0:5 lm. This surface is technically smooth, hence the effect of surface roughness on frictional losses is assumed to be negligible. The procedure for evaluating the heat transfer coefficients is described in [5]. The total relative uncertainties for pressure drop are below 2% and for the heat transfer coefficients below 3%, at a confidence level of 95% (cf. [4]). Systematic uncertainties in temperature, pressure and flow measurements were minimized by calibration of the measuring chains, covering the whole range of operation. Fig. 4. Flow sheet of the experimental facility used for the measurement of pressure loss and evaluation of heat transfer coefficient inside pillow plates (see also [5]).

5. Experimental setup The experimental facility used to measure pressure loss and determine heat transfer coefficients in pillow plates is shown in Fig. 4. Its detailed description is given in [4,5]. The pressure loss used for the validation of the CFD simulations was measured between two measuring ports (ports 1 and 2 in Fig. 4) placed far enough from the inlet and outlet. In this way, only the fully developed flow region was considered. The ports were made by drilling bore holes of 2 mm in diameter through the pillows, directly at the intersection of the longitudinal distance and transversal distance between the welding spots. The distance

6. Validation In Fig. 5(a), specific pressure loss determined by CFD simulations is compared with experimental values obtained by our group. The numerical and experimental results agree well over the entire studied Reynolds numbers range (1000 6 Re 6 8000). The relative deviation is below 5%, with the exception of Re ¼ 1000, at which it is 8%. It can be concluded that the complex steady-state fluid dynamics in the pillow-plate channel is well captured by the CFD model. For the comparison of Nusselt numbers, we used a constant wall heat flux boundary condition (q_ w ¼ const) instead of the constant wall temperature mentioned above. This was done, because in the experiments the pillow-plate wall was heated electrically (cf. [4]). Nusselt numbers determined using

Fig. 5. Comparison of numerical and experimental results for (a) turbulent pressure loss, (b) Nusselt numbers for Pr ¼ 6 and Re ¼ 1000; 2000; 4000; 6000; 8000, (c) Nusselt numbers for Re ¼ 4000 and Pr ¼ 5; 20; 50; 100.

521

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

q_ w ¼ const are approximately 6% higher than those determined with T w ¼ const. For fully developed turbulent flow in pipes and other simple channel geometries, no difference in Nusselt numbers determined either with constant heat flux or constant temperature at the walls can be found (see, e.g., [18]); in pillow plates the deviation possibly originates from the recirculation zones, which are heated up more strongly when q_ w ¼ const applies. This is illustrated later in Fig. 7(a). In Fig. 5(b) and (c), a comparison of simulated and experimentally determined Nusselt numbers is shown. Fig. 5(b) gives it for Pr ¼ 6 and varied Reynolds number, whereas Fig. 5(c) for Re ¼ 4000 and varied Prandtl number. In both figures, a good agreement is visible. At the lowest Reynolds number of this study, Re ¼ 1000, the flow is still fully turbulent. At the largest Reynolds number, Re ¼ 8000, the specific pressure loss exceeds 1 bar=m, which results in an inefficient operation of the heat exchanger. 7. Results and discussion 7.1. Fluid dynamics and heat transfer in pillow plates The presence of welding spots in pillow-plate channels as well as the waviness of the channel walls lead to a strong deflection of the flow and to the rise of pronounced secondary-flow effects, as shown in Fig. 6. In the immediate vicinity of these spots, the channel is narrowest, leading to a large flow resistance. As the result, the fluid is guided away from the welding spots in radial direction, thus producing large regions of recirculating fluid behind them. The primary flow is strongly deflected by the welding spots and the reciruculation zones and follows a nearly sinusoidal path. In Fig. 7(a) the steady vortices in the wake of the welding spots are shown in more detail. They have a weak ‘‘tornado-type” character. Fluid near the wall is captured by the vortex core and transported normally to the wall, up to the symmetry plane (z ¼ 0). The unit vectors of the wall shear stress (s=jsj) in Fig. 7(a) and the streamlines adjacent to the wall (blue) elucidate the inward fluid flow (near the wall) towards the recirculation zones. In the symmetry plane, the streamlines at the boundary of the recirculation zones follow a slightly different path than those close to the wall. Consequently, the rise of weak cork-screw vortices at the boundaries of the recirculation zones can be observed. These vortices, however, do not seem to cause any significant increase in wall shear stress or heat transfer rate. Fig. 7(b) shows, based on the scalar velocity field, the characteristic ‘‘two-zone” flow pattern encountered in pillow plates. Zone 1

Fig. 6. Illustration of the characteristic flow pattern in pillow plates represented by streamlines. The bottom wall of the channel is solid while the upper wall is transparent for the sake of visualization.

is represented by the meandering core flow and zone 2 by the recirculation zones. The size and shape of these zones depends on the pillow-plate geometry (e.g. welding spot pattern and inflation height) and will be discussed further. As shown in [3], flow separation is caused by an adverse pressure gradient; it appears downstream of the stagnation point of the welding spot at approximately two-thirds of the spot circumference. Furthermore, in Fig. 7(b), large velocity gradients are observed at the boundary between the meandering core and the recirculation zones, where the local velocity increases abruptly by a factor of about 11. The influence of the flow phenomena shown in Fig. 7(b) on heat transfer is reflected by the field of the normalized wall heat flux in Fig. 7(c). Heat transfer in the recirculation zones is clearly poor; this results in hot spots in the temperature field (see Fig. 7(d)), especially at the eddy centers and in the immediate vicinity of the welding spots. Consequently, heat is transferred mainly in the meandering core region, where the local Reynolds number is approximately twice the mean Reynolds number of the total channel flow. The periodic acceleration of the flow in stream-wise direction of the meandering core (cf. Fig. 7(b)), leads to a local increase of the wall heat flux (cf. Fig. 7(c)), caused by local reduction of the boundary layer thickness. The simulations showed that the mean heat transfer coefficient in the pillow-plate channel with the geometry presented in Fig. 7 was proportional to Re0:79 . The exponent of the mean Reynolds number is mainly determined by the meandering core and is close to the value 0.8, which is typical for turbulent heat transfer in smooth pipes (see, e.g., [19]). Consequently, from the thermohydraulic point of consideration, the meandering region resembles a typical channel flow, which is bounded by the pillow-plate walls and the recirculation zones. Furthermore, the simulations revealed that the heat transfer coefficient for the geometry shown in Fig. 7 is proportional to Pr0:38 , which is close to the exponent 1=3 typical for turbulent heat transfer in channels. Steimle [20] analyzed the exponents of the Reynolds and Prandtl numbers in correlations for turbulent forced convection heat transfer in numerous geometries and found that in most cases, the exponent of the Reynolds number is twice that of n

the Prandtl number, i.e. Nu ¼ CðRe2 PrÞ . The factor C and the exponent n depend on the geometry. Considering the exponents for Re and Pr obtained from the simulations in this study, the expression of Steimle [20] also holds for pillow plates. 7.2. Influence of pillow-plate geometry on flow characteristics The ‘‘two-zone” character of the flow pattern shown in Fig. 7 is characteristic for the flow in pillow-plates; however, depending on the pillow-plate geometry, the size and shape of the recirculation zones can vary. In the simulations we observed three subcategories of the characteristic flow pattern, namely for pillow plates with sr > 1:56 (longitudinal type), sr 6 1 (transversal type) and 1 < sr 6 1:56 (mixed type). Both the two-zone character of the flow and the differences between the flow patterns of these categories can be quantified with the help of the following parameters:

wR ¼

fR fDp

ð11Þ

wA ¼

ARZ Atot

ð12Þ

wQ ¼

Q_ RZ Q_ tot

ð13Þ

wRe ¼

Remc Retot

ð14Þ

522

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

Fig. 7. Vortex structures in the wake of the welding spots (a), developed velocity field (b), field of normalized wall heat flux (c), and temperature field (d) for turbulent flow in a pillow-plate channel for Re ¼ 2000. Geometry: sT ¼ 42 mm, sL ¼ 36 mm, dsp ¼ 10 mm and di ¼ 6 mm.

The parameter wR expresses the portion of the pressure loss coefficient fDp responsible only for skin friction fR . The size of the recirculation zones is characterized by wA , which is defined as the ratio of the wall area covered by the recirculation zones to the total area of the wall. The thermal activity of the recirculation zones is quantified by wQ , and is defined as the ratio of the rate of heat transferred only in the recirculation zones to the total rate of heat transferred normal to the wall. 7.2.1. Heat and fluid flow characteristics of longitudinal-type pillow plates (sr P 1:56) The results presented in Fig. 7 are for longitudinal-type pillow plates. They are characterized by recirculation zones, which have a ‘‘flame-like” shape and have a length close to sL . Due to this shape, the primary flow is deflected significantly and follows a ‘‘meandering-like” path. The profile of wA in Fig. 8(a) shows that the flow pattern in pillow-plates with sr > 1:56 is only weakly dependent on the Reynolds number. Thus, the area of the recirculation zones varies only slightly with Re and occupies approximately 45% of the total wall area. These regions cause the major portion of the form drag, which is responsible for more than 50% of pressure loss in the pillow-plate channel, as can be confirmed by the wR profile in Fig. 8(a). The larger the Reynolds number is, the larger is the size of the recirculation zones and, hence, the contribution of form drag to the total pressure loss. Note that the form drag coefficient in longitudinal-type pillow plates is almost constant for all Reynolds numbers (Fig. 8(b)), so that the trend of the pressure loss

coefficient is similar to the trend of the fanning friction factor (cf. Fig. 8(b)). The rate of heat transferred normal to the wall in the recirculation zones is only 15% of the total, even at the highest Reynolds number, although almost 50% of the wall area is covered by recirculation zones. This confirms the poor thermal activity in the recirculation zones and that, conversely, the meandering core dominates heat transfer. The latter is mainly related to the large Reynolds number in this region (wRe  2). 7.2.2. Heat and fluid flow characteristics of transversal-type pillow plates (sr 6 1) The flow pattern shown in Fig. 9 is characteristic for transversal-type pillow plates with sr 6 1. In contrast to longitudinal-type, the recirculation zones in transversal-type pillow plates extend to the next downstream welding spot and thus occupy a larger portion of the channel. The shape of the recirculation zones affect the flow path of the meandering core; the latter is deflected less significantly by the welding spots and follows a more linear path. The flow pattern in transversal-type pillow plates is also largely independent of the Reynolds number, as can be seen from the wA curve in Fig. 10(a). In contrast to longitudinal-type pillow plates, the values of wA in Fig. 10(a) show that the recirculation zones in transversal-type pillow plates are larger. However, the contribution of the form drag coefficient to the overall pressure loss coefficient is smaller than for

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

523

Fig. 8. Dimensionless two-zone variables (a) and friction coefficients (b) as functions of the Reynolds number for Pr ¼ 6. Geometry: sT ¼ 42 mm, sL ¼ 36 mm, dsp ¼ 10 mm and di ¼ 3 mm.

Fig. 9. Vortex structures in the wake of the welding spots (a) and developed velocity field (b) for turbulent flow in transversal-type pillow plates for Re ¼ 2000. Geometry: sT ¼ 72 mm, sL ¼ 21 mm, dsp ¼ 10 mm and di ¼ 6 mm.

Fig. 10. Dimensionless two-zone variables (a) and friction coefficients (b) as functions of the Reynolds number for Pr ¼ 6. Geometry: sT ¼ 72 mm, sL ¼ 21 mm, dsp ¼ 10 mm and di ¼ 3 mm.

longitudinal-type pillow plates and, moreover, smaller than the contribution of the fanning friction factor (wR > 0:5; cf. Fig. 10(b)). This is related to the shape of the recirculation zones, which occupy less channel cross-section. As the consequence, the meandering core is less deflected and, hence, the resistance to flow is smaller. The thermal activity of the recirculation zones in transversaltype pillow plates is poorer than in longitudinal type. Whilst the recirculation zones in the former type are larger in size, the portion of heat transferred in these zones (wQ ) is similar to that in longitudinal type. The mean Reynolds number in the meandering core is larger than in longitudinal-type pillow plates (cf. Fig. 10(a)) because of the more homogeneous velocity profile, related to the larger channel cross-section.

7.2.3. Heat and fluid flow characteristics of mixed-type pillow plates (1 < sr 6 1:56) The flow patterns in longitudinal-type and transversal-type pillow plates are largely independent of the Reynolds number, i.e. wA varies only slightly with Re. In mixed-type pillow plates, however, the flow pattern varies substantially with Re. As can be seen in Fig. 11, at small Reynolds numbers (Re 6 2000) the flow pattern resembles that of longitudinal-type pillow plates. At larger Reynolds numbers (Re > 2000), the recirculation zone increases substantially in size until it reaches the next downstream welding spot. At Re ¼ 8000, the flow pattern in Fig. 11 combines features of both longitudinal-type and transversal-type pillow plates. The variation of the flow pattern with the Reynolds number can be followed by the trend of wA in Fig. 12(a). Up to Re ¼ 4000; wA

524

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

increases strongly and thereafter only weakly with the Reynolds number. A similar behavior of wQ is visible, since this parameter is directly related to the size of the recirculation zones. The dependency between the size of the recirculation zones, i.e. wA , and form drag is shown in Fig. 12(b). Due to the strong increase in the size of the recirculation zones for Re < 4000, a large increase in the form drag coefficient fD is observed. The variation of the pressure loss coefficient fDp with the Reynolds number shown in Fig. 12(b) is characteristic for mixed-type pillow plates. The profile of fDp in Fig. 12(b) can be explained by the sum: fDp ¼ fR þ fD . The aforementioned increase of the form drag coefficient at small Reynolds numbers causes a sudden increase in fDp , which then gradually becomes constant towards larger Re values. 7.3. Thermo-hydraulic efficiency of pillow-plates The thermo-hydraulic efficiency is represented by the ratio of the benefit (the heat transferred) and the effort (pumping power) required:

Q_

hA DT

e ¼ _ ¼ _w W V Dp

ð15Þ

The choice of characteristic geometrical parameters of pillow plates, e.g. the welding spot pattern and the channel height, strongly influences their thermo-hydraulic efficiency. When

Fig. 11. Developed velocity field of turbulent flow in mixed-type pillow plates for different Reynolds numbers (at symmetry plane ðx; y; 0Þ). Geometry: sT ¼ 42 mm, sL ¼ 30 mm, dsp ¼ 10 mm and di ¼ 3 mm. The red cross marks the end of recirculation zones.

designing pillow plates, it is important to quantify this efficiency in order to (a) obtain the optimal geometry and (b) guarantee the most economical operation of the heat exchanger for a given task. Eq. (15) is used to quantify the efficiency of the pillow-plate geometries investigated in this work. For a more convenient comparison with other geometries, the efficiency e is normalized:

e ¼

e eref

ð16Þ

so that it ranges from 0 to 1. Normalization is performed by dividing Eq. (15) by the efficiency of the geometry 42/72/6/10 (eref ) at Re ¼ 1000. The latter showed the highest efficiency of all pillow plates investigated in this study excluding the optimization study which will be presented later. Fig. 13(a) shows an example of the variation of e with the Reynolds number for the geometry 42/72/6/10. All other pillow plate geometries demonstrated a similar trend. The largest efficiency is obtained at the lowest Reynolds number Re ¼ 1000; it decreases asymptotically thereafter. The latter is generally due to the different dependencies of the heat transfer coefficient and pressure loss on the Reynolds number. For example, the pillow plate 42/72/6/10 showed that h Re0:76 and that

Dp Re1:63 . Hence, compared to the heat transfer coefficient, pressure loss increases more significantly with the Reynolds number and thus decreases the efficiency. Fig. 13(b) shows how the welding spot pattern and the channel height influence the efficiency. From left to right, i.e. from transversal-type towards longitudinal-type pillow plates, a minimum in the efficiency, attributed to the equidistant pattern sr ¼ 1, is observed. Hence, equidistant pillow plates produced the smallest efficiency of all pillow plates studied. Considering that each welding spot represents a resistance to flow, in the case of sr ¼ 1, the number of welding spots per square meter is largest and so is the resistance to flow. In the case of longitudinal-type pillow plates, the efficiency is better than for sr ¼ 1 (for the same sT ), because the number of welding spot rows in stream-wise direction is smaller. For sr  1, the efficiency e seems to approach a constant value. Hence, welding spot patterns, whose ratio sr exceeds that of the typical triangular arrangement with sr > 1:9, do not seem to provide further advantages. On the other hand, approaching more ‘‘transversal” geometry (i.e. with decreasing sr ) yields better efficiency. This is because the portion of form drag (‘‘bad pressure loss”) to overall pressure loss is reduced compared to longitudinal-type pillow plates. Form drag increases overall pressure loss more significantly without contributing to heat transfer. The geometry 42/42/6/10 is a good example. The heat transfer coefficient shows that h Re0:74 and

Fig. 12. Dimensionless two-zone variables (a) and friction coefficients (b) as functions of the Reynolds number for Pr ¼ 6. Geometry: sT ¼ 42 mm, sL ¼ 30 mm, dsp ¼ 10 mm and di ¼ 3 mm.

525

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

Fig. 13. Thermo-hydraulic efficiency of (a) the transversal type pillow-plate 42/72/6/10 as a function of the Reynolds number, and (b) as a function of sr and inflation height for Re ¼ 2000.

that Dp Re1:9 . It becomes apparent, that the difference in the exponents of the Reynolds number for h and Dp is larger than for e.g. the transveral type pillow plate 42/72/6/10. Thus, the increase in Dp with Re, compared to the increase of h with Re, is stronger than for pillow plates with smaller sr . For sr ! 0, the efficiency should asymptotically approach the value of a plane channel. Furthermore, the efficiency increases with the channel height (see Fig. 13(b)). The latter increases the mean channel crosssection, and, consequently, the mean stream velocity decreases for constant Reynolds number. As the result, pressure loss is smaller. In summary, the optimal design can be obtained with a pillow plate having a small sr , a large channel height and operating it at ‘‘lower” Reynolds numbers. Although the smallest Reynolds number produces the highest efficiency, it delivers the smallest heat transfer coefficient. This means, that the heat exchanger will become large in size (large heat transfer area) to compensate for the small h, thus resulting in high equipment cost. A reasonable compromise between thermo-hydraulic efficiency and heat transfer coefficient is achieved by using the Reynolds number at the intersection between the two tangents, as illustrated in Fig. 13(a).

8. Geometry optimization As discussed above, the recirculation zones observed in the wake of the welding spots increase pressure loss and are less effective for heat transfer. Consequently, by reducing the size of these zones, an improvement of the thermo-hydraulic performance of pillow plates can be achieved. For this reason, we investigated the effect of oval-shaped welding spots, which are more streamlined than round ones, on the recirculation zones.

Fig. 14 shows an example of two pillow-plate geometries, which have the same welding spot pattern and inflation height, but different welding spot shapes. In this example, the transversal welding spot diameters dSP of both geometries are equal; however, the longitudinal diameter of the oval welding spots is two times the transversal diameter. This changes both the shape of the welding spots and the waviness of the channel walls (cf. Fig. 14); hence, also the local cross-sectional areas of the channel are changed. Table 2 summarizes the pillow-plate geometries investigated in this optimization study. The letter ‘‘D” at the end of the abbreviation in Table 2 denotes pillow plates with oval welding spots, which have the same transversal diameter as that of round ones. Accordingly, the letter ‘‘A” represents pillow plates with oval welding spots, which have the same surface area as that of round spots. The effect of the welding spot shape on the size of the recirculation zones in longitudinal type pillow plates is presented in Fig. 15. The use of oval welding spots clearly reduces the size of the secondary flow (cf. Fig. 15(b) and (c)). Although the shape of these welding spots is more streamlined, they do not lead to a shift of the separation point. Rather the reduction in the size of the recirculation zones is caused by the channel waviness. The latter is notably different compared to a pillow plate with round welding spots. The width of the diagonal cross-section of the channel with elliptical welding spots (bD;E ) is smaller than for the pillow plate with round spots (bD;E < bD;C ). This leads to an earlier redirection of the meandering core, due to the opposing action of the next (diagonally placed) downstream welding spot on the flow. Consequently, the recirculation regions are ‘‘squashed” and thus reduced in size. This results in the meandering core following a more linear path with less deflection, compared to behavior shown in Fig. 15(a). In Fig. 16, the simulation results for the optimization of longitudinal-type pillow plates are summarized. As can be seen, the geometry 42/72/6/D brings a reduction in wA by 23%, which indicates a significant decrease in the size of the recirculation

Table 2 List of investigated pillow-plate geometries for the optimization study.

Fig. 14. Pillow plate with round welding spots (a) and with elliptical welding spots (b).

2sL (mm)

sT (mm)

di (mm)

dSP (mm)

lSP =dSP (mm)

Comments (abbreviation)

Abbreviation

72 72 72 72 42 42 42 42

42 42 42 42 72 72 72 72

6 6 6 6 6 6 6 6

10 10 7.07 7.2 10 10 7.07 7.2

1 2 2 1 1 2 2 1

– – Around ¼ Aov al – – – Around ¼ Aov al –

72/42/6/10 72/42/6/D 72/42/6/A 72/42/6/7.2 42/72/6/10 42/72/6/D 42/72/6/A 42/72/6/7.2

526

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

Fig. 18. Results of the optimization study. Thermo-hydraulic characteristics of transversal type pillow-plates with different welding spot geometries, compared to the reference pillow plate 42/72/6/10. The results are based on the same volumetric flowrate in all geometries.

Fig. 15. Developed velocity field in a pillow plate with the geometry 72/42/6/10 (a), 72/42/6/D (b), 72/42/6/A (c) and 72/42/6/7.2 (d). The simulations (a)–(d) were performed for the same volumetric flowrate. For geometry (a), this flowrate corresponds to Re ¼ 2000.

Fig. 16. Results of the geometry optimization study. Thermo-hydraulic characteristics of longitudinal type pillow-plates with different welding spot geometries, compared to the reference pillow plate 72/42/6/10. The results are based on the same volumetric flowrate in all geometries.

zones. As a result, the portion of form drag in total resistance is also reduced (cf. increase in wR ), which leads to a decrease in pressure loss of up to 12%. On the other hand, the product hAw is reduced by almost the same amount, thus reducing the improvement of the thermo-hydraulic efficiency e. The product hAw is more informative than just h, since the heat transfer area in 42/72/6/D is reduced by 5%, due to the larger surface of the oval welding spots compared to round ones. Therefore, we also investigated oval welding spots, which have the same surface area as the round ones in 42/72/6/10. The flow pattern in such pillow plates is seen in Fig. 15(c), and is similar to that in Fig. 15(b). The geometry

42/72/6/A results in a significant reduction in pressure loss of up to 35%, while hAw decreases by only 9%. The combined result is an improvement in efficiency e by as much as 37%. Since the transversal diameter of the smaller oval welding spot in 42/72/6/ A leads to a smaller form drag, we also investigated the use of round welding spots with a smaller diameter than in 42/72/6/10. The flow pattern of this case is shown in Fig. 15(d). The small diameter round welding spots in 42/72/6/7.2 are even more efficient than the small oval welding spots (42/72/6/A). They lead to an efficiency improvement by as much as 43% (cf. Fig. 16). Since the use of oval welding spots and small round ones improves the overall thermo-hydraulic performance of longitudinal-type pillow plates, we studied their impact on transversal-type pillow plates, too. The results of this study are presented in Figs. 17 and 18. In contrast to longitudinal-type pillow plates, the flow pattern in transversal-type pillow plates with oval welding spots does not change notably. The size of the reciruculation zones is reduced (cf. Fig. 17(b)); this time mainly due to the smaller distance lR between successive welding spots. Interestingly, the efficiency in the transversal-type geometry 42/72/6/D is even reduced by 5%, compared to the reference geometry 42/72/6/10. On the other hand, by using oval welding spots with the same surface area as the round ones in 42/72/6/10, i.e 42/72/6/A, the thermo-hydraulic performance is improved by 22%. Similar to longitudinal-type pillow plates, the best improvement in transversal-type pillow plates is achieved by using small round welding spots (42/72/6/7.2), with an efficiency of up to 25% (cf. Fig. 18). In summary, the use of oval welding spots or small round ones can significantly improve the thermo-hydraulic efficiency of pillow plates, especially of the longitudinal type.

Fig. 17. Developed velocity field in a pillow plate with the geometry 42/72/6/10 (a), 42/72/6/D (b), 42/72/6/A (c) and 42/72/6/7.2 (d). The simulations (a)–(d) were performed for the same volumetric flowrate. For geometry (a), this flowrate corresponds to Re ¼ 2000.

M. Piper et al. / International Journal of Heat and Mass Transfer 94 (2016) 516–527

527

9. Conclusions

References

In this study we performed a comprehensive CFD analysis of fluid dynamics and heat transfer in pillow plates. The simulations showed that the fully developed turbulent flow in pillow plates is characterized by two distinctive regions, namely, (a) recirculation zones arising in the wake of the welding spots and, (b) the core flow, which is bounded by the walls and the recirculation zones. The latter flow represents a primary effect while the former is secondary. The secondary flow produces form drag, which increases pressure loss in the channel. Depending on the geometry, we identified three sub-categories of the characteristic flow pattern, a longitudinaltype, transversal-type and mixed-type pillow plates. They differ mainly in the size and shape of the recirculation zones. Whereas the characteristic flow pattern in longitudinal-type and transversal-type pillow plates is only weakly dependent on the Reynolds number, in mixed-type pillow plates, the recirculation zones increase significantly in size with increasing Reynolds number. We also investigated the thermo-hydraulic efficiency of pillow plates. The largest efficiencies were observed at the lowest Reynolds number and decreased rapidly at larger Reynolds numbers. This effect is caused by the fact that, with increasing Reynolds number, pressure loss grows more rapidly than the heat transfer coefficient. Furthermore, the lowest efficiencies were observed for pillow plates, in which the longitudinal and transversal welding spot pitches were equal. The highest efficiencies were obtained for transversal-type pillow plates with the larger inflation height. Finally, we performed a geometry based optimization study for improving the thermo-hydraulic efficiency of pillow plates. We found that by using ‘‘oval-shaped” welding spots, which are more streamlined than round ones, the efficiency can be increased by as much as 37% for longitudinal-type pillow plates and by 22% for transversal-type pillow plates. Greater improvements in the efficiency were found for reduced size of the welding spots. The results of this study offer detailed information about heat and fluid flow characteristics in pillow plates, which are useful for the development of more detailed and reliable design methods. The characteristic flow pattern identified here can be used for the development of ‘‘flow-pattern oriented” design methods for pillow plates, e.g. a ‘‘2-zone-model”, in which the meandering core and the recirculation zones are modeled separately.

[1] J. Mitrovic, R. Peterson, Vapor condensation heat transfer in a thermoplate heat exchanger, Chem. Eng. Technol. 30 (7) (2007) 907–919, http://dx.doi.org/ 10.1002/ceat.200700082. [2] J. Mitrovic, B. Maletic, Numerical simulation of fluid flow and heat transfer in thermoplates, Chem. Eng. Technol. 34 (9) (2011) 1439–1448, http://dx.doi.org/ 10.1002/ceat.201100271. [3] M. Piper, A. Zibart, J.M. Tran, E.Y. Kenig, A numerical study on turbulent singlephase flow and heat transfer in pillow-plates, in: Proceedings of the 15th International Heat Transfer Conference, Kyoto, Japan, 2014. http://dx.doi.org/ 10.1615/IHTC15.hex.008929. [4] J.M. Tran, S. Sommerfeld, M. Piper, E.Y. Kenig, Investigation of pillow-plate condensers for the application in distillation columns, Chem. Eng. Res. Des. 99 (2015) 67–74, http://dx.doi.org/10.1016/j.cherd.2015.03.031. [5] J.M. Tran, M. Piper, E.Y. Kenig, Experimentelle Untersuchung des konvektiven Wärmeübergangs und Druckverlustes in einphasig durchströmten Thermoblechen, Chem. Ing. Tech. 87 (3) (2015) 226–237, http://dx.doi.org/ 10.1002/cite.201400140. [6] B. Maletic, Modelling and Numerical Analysis of Fluid Flow and Heat Transfer in Thermoplates, University of Paderborn, 2009 (Ph.D. thesis). [7] M. Piper, A. Olenberg, J.M. Tran, E.Y. Kenig, Determination of the geometric design parameters of pillow-plate heat exchangers, Appl. Therm. Eng. 91 (2015) 1168–1175, http://dx.doi.org/10.1016/j.applthermaleng.2015. 08.097. [8] T.H. Shih, W.W. Liou, A. Shabbir, Z. Yang, J. Zhu, A new k— eddy viscosity model for high reynolds number turbulent flows, Comput. Fluids 24 (3) (1995) 227–238, http://dx.doi.org/10.1016/0045-7930(94)00032-T. [9] P.A. Durbin, B.A.P. Reif, Statistical Theory and Modeling for Turbulent Flows, John Wiley and Sons, 2003. [10] A. Lautenschleger, A. Olenberg, E.Y. Kenig, A systematic CFD-based method to investigate and optimise novel structured packings, Chem. Eng. Sci. 122 (2015) 452–464, http://dx.doi.org/10.1016/j.ces.2014.09.040. [11] S.A. Owens, M.R. Perkins, R.B. Eldridge, Computational fluid dynamics simulation of structured packings, Ind. Eng. Chem. Res. 52 (2013) 2032– 2045, http://dx.doi.org/10.1021/ie3016889. [12] S. Jain, A. Joshi, P.K. Bansal, A new approach to numerical simulation of small sized plate heat exchangers with chevron plates, J. Heat Transfer 129 (2007) 291–297, http://dx.doi.org/10.1115/1.2430722. [13] Q. Wang, Q. Chen, G. Chen, M. Zeng, Numerical investigation on combined multiple shell-pass shell-and-tube heat exchanger with continuous helical baffles, Int. J. Heat Mass Transfer 52 (2009) 1214–1222, http://dx.doi.org/ 10.1016/j.ijheatmasstransfer.2008.09.009. [14] W. Rodi, Experience with two-layer models combining the k— model with one-equation model near the wall, AIAA (1991) 91–0216. [15] F. Lien, W.L. Chen, M.A. Leschziner, Low-Reynolds-number eddy-viscosity modelling based on non-linear stress-strain/vorticity relations, in: Proc. 3rd Symp. on Engineering Turbulence Modelling and Measurements, 1996, pp. 91–100. [16] M. Wolfshtein, The velocity and temperature distribution in one-dimensional flow with turbulence augmentation and pressure gradient, Int. J. Heat Mass Transfer 12 (3) (1969) 301–318, http://dx.doi.org/10.1016/0017-9310(69) 90012-X. [17] P.A. Davidson, Turbulence – An Introduction for Scientists and Engineers, Oxford University Press, New York, 2004. [18] W. Kays, M. Crawford, B. Weigand, Convective Heat and Mass Transfer, McGraw-Hill, 2005. [19] H. Kraussold, Die Wärmeübertragung an Flüssigkeiten in Rohren bei turbulenter Strömung, Forsch. Ingenieurwes. 4 (1) (1933) 39–44. [20] F. Steimle, Der allgemeine Zusammenhang zwischen Wärmeübergang und Druckabfall turbulenter Strömungen und seine Auswirkung auf die Formulierung von Funktionen der Nusselt-Zahl, Verlag C.F. Müller Karlsruhe, 1970.

Acknowledgments The authors sincerely acknowledge the financial support from the German Federal Ministry of Education and Research within the scope of the joint project InnovA2 (Grant No. 033RC1013E), and from the EC FP7 project ENER/FP7/296003165 ‘‘Effcient Energy Integrated Solutions for Manufacturing Industries” (EFENIS).