Guy Clarence Semassou

Abdel Deen

Derrick Vital

Dai Tometin
email: derricktometin@gmail.com

Roger Ahouansou

• Kouamy

Victorin Chegnimonhan
email: victorin.chegnimonhan@gmail.com

Tognon Clotide Guidi
email: guidi65@mail.ru

* Guy

Clarence Semassou
email: clarence.semassou@gmail.com

Contribution to the optimization of the autonomous photovoltaic solar systems with hybrid storage for loads with peak power with constraints of volume and loss of power supply probability

Keywords: Optimization, Lifetime cost, PV solar system, Ultracapacitor C Cap Nominal ultracapacitor capacitance

published or not. The documents may come

Contribution to the optimization of the autonomous photovoltaic solar systems with hybrid storage for loads with peak power with constraints of volume and loss of power supply probability

Project's lifetime

Introduction

The importance attached to optimization is increased for the autonomous generators of renewable energy, because the safety of a connection to the national network for the backup supply does not exist. A certain number of diferent methods are used for optimization of the intermittent autonomous generators. These methods include intuitive analytical methods and approaches of simulation of system [START_REF] Semassou | AIDE A LA DECISION POUR LE CHOIX DE SITES ET SYSTEMES ENERGETIQUES ADAPTES AUX BESOINS DU BENIN[END_REF].

A dimensioning method of a PV generator is described in the guide of the IEEE for dimensioning of panels and batteries in the autonomous photovoltaic systems [START_REF]IEEE Guide for Array and Battery Sizing in Stand-Alone Photovoltaic (PV) Systems[END_REF]. The number of required photovoltaic panels is based on the solar radiation, the ratio network/load, the losses of the system and the load. If the load is constant, the solar radiation of the most unfavorable month is used in calculations. If the load varies, the most unfavorable case is selected like month with weakest ratio network/load. The number of necessary batteries is given according to the days of autonomy (DOA), by supposing that the photovoltaic generator does not provide any power.

Lagorse et al. [START_REF] Lagorse | Sizing optimization of a stand-Alone street lighting system powered by a hybrid system using fuel Cell, PV and battery[END_REF] developed a method of optimization for a PV system of lighting with storage of hydrogen and battery. Methodology suggested uses an approach in two stages to determine the optimal solution. Initially a genetic algorithm is used to ind a tendency of convergence; the value obtained starting from the genetic algorithm is used like starting point for an algorithm simplex. The objective is to minimize the cost on cycle of life of the system. The parameters taken into account during optimization are the power of the PV and the fuel cell, the angle of inclination of the PV system, the capacity of the battery, the fuel cell of starting and stop and the state of charge SOC of the battery. The current and the tension of the battery were not taken into account during the process.

Ekren and Ekren [START_REF] Ekren | Simulation based size optimization of a PV/wind hybrid energy conversion system with battery storage under various load and auxiliary energy conditions[END_REF] carried out an optimization based on the simulation of a system wind-PV with storage of battery for a mobile station of communication GSM and for a domestic proile. Optimization was carried out with various proiles of load and auxiliary energy costs. Simulations are carried out with software ARENA and Opt-Quest. Opt-Quest requires a maximum and minimum number of components to enter for optimization to carry out.

Yang et al. [START_REF] Yang | Optimal sizing method for stand-Alone hybrid solar-Wind system with LPSP technology by using genetic algorithm[END_REF] developed a methodology of dimensioning of an autonomous hybrid system wind-PV with storage of battery, based on the genetic algorithms. Optimization is carried out by taking account of the loss of probability of power supply (LPSP) and of the annualized cost of the system. Methodology describes the number of PV panels, batteries and wind mills as well as the angle of inclination of the panels and the height of the wind mills. The genetic algorithm uses time data and models of components to seek the optimal solution. The initial values for each component are necessary for the process of optimization.

Borowy and Salameh [START_REF] Borowy | Methodology for optimally sizing the combination of a battery bank and PV array in a wind/ PV hybrid system[END_REF] developed a methodology to dimension in an optimal way the combination of a park of batteries and a PV generator in a hybrid system wind-PV autonomous. Methodology used makes it possible to minimize the economic cost of the system while ensuring without interruption the cover of the requirements in electric power. In methodology, the size of the wind mill is deined, with an iterative approach used to determine the number of batteries and necessary PV panels according to constraint LPSP. Optimization was carried out for each 1 h of typical day of each month. Seeling-Hochmuth [START_REF] Seeling-Hochmuth | A combined optimisation concept for the design and operation strategy of hybrid-PV energy systems[END_REF] jointly optimized the strategy of size and control for the hybrid systems, wind-PV, coupled with batteries with an auxiliary diesel generator. Optimization is carried out in two stages, the principal algorithm treats dimensioning of the system while the second subsystem treats optimization of control. For the size of selected system, the limits of control are adjusted to minimize the objective function subjected to constraints of exploitation, for example to minimize the cost of the system. The result of the process is an optimized coniguration system and a strategy of control for the application and the site given.

The methods of optimization discussed in the literature make the average of the powers generated and consumption on an hourly basis to dimension the components of the solar system. The use of the average time data is satisfactory when in fact the requirements in energy are taken into account. However, the amplitude of the peaks of power is attenuated when it is realized on an hourly basis. Consequently, when we consider the efects of the peaks power on the system of storage, shorter temporal periods are necessary to collect the efects of the peaks power on the tension of the battery. Optimization by taking account peaks of power constitutes the principal objective of our research. In addition to that, we propose a method to correctly size the combination of the battery and the supercondensator in order to avoid the excessive costs of the system while having a constraint on the volume of the system.

Materials and methods

Adopted proiles of consumption

A certain number of proiles of load were developed to analyze the method of optimization under variable conditions. The proiles of load used were developed to represent a peak, an impulse and a domestic proile of load. The proile of peak load was developed to represent the starting of an engine with D.C. current. The current necessary to the initial starting of the engine is deinitely higher than the current necessary to maintain the engine in steady operation. Typically, the inrush current of the engine could be six times the current in steady operation [START_REF] Kaplanis | The Design and Integration of Possible PV Conigurations to Determine the Most Cost Efective System for a Household[END_REF].The proile of peak load has a power peak of 200 W and the power continues is of approximately 32 W.The engine is used during 45 min every hour, as illustrated on Fig. 1a, the proile of load being repeated daily over the period of optimization.

A proile of load of impulse is the second developed to evaluate the method. The load of impulse functions during 200 s over 250 s with a power of high impulse of 50 W and a basic power lower by 10 W the cyclic report/ratio of the impulses is 0.5 with one 20 s period. The proile of load is illustrated on Fig. 1b and functions uninterrupted.

The domestic proile of load was used to examine the method of optimization under a variable proile of load. The proile of domestic load annual describes on Fig. 2 is considered for our study. The proile has variable peaks of power occurring each day.

Variation of the sunning and the temperature according to time

The program of optimization developed requires the solar radiation on a horizontal level and the ambient temperature like entries to generate the power of maximum exit of PV panels over the duration of optimization. As the solar radiation varies during the year, a proile of annual solar radiation was obtained from a point located at Benin, precisely in

Cotonou. The data weather of the ASECNA were considered. Figure 3a, b give respectively a chart of the solar radiation and ambient temperature recorded on the site.

Description of the energy system

The studied system includes various elements whose photovoltaic solar panels, the controller of load MPPT, the system of storage of energy (battery VRLA), the supercondensator, the controller of load for hybrid storage, and load (Fig. 4).

Modeling

Model of determination of the numbers minimal and maximum of panels

In autonomous PV systems, PV panels are the only source of generation of energy. Consequently, PV panels must be able to generate the energy minimum required by the power consumption of the loads by taking account of the energy losses of PV system. As the power produced by the photovoltaic generator varies throughout the year, it exists a compromise between the increase in the number of photovoltaic panels and the increase in the system of storage, which feeds the load during the periods of weak or any solar radiation. To ind the quantity optimal of panels photovoltaic to feed the load throughout the year with a System with Hybrid Storage, the number of photovoltaic panels varies between minimal and maximum. The minimal number of photovoltaic panels required to fulill the where the day varies from 1 to 365, so that the program analyzes each day of the year

Model of capacity of the supercondensator

In the System with Hybrid Storage, the park of supercondensator is used to provide the needs for power of point to the load which discharges the system from storage with a low average power provided by the bench of batteries. The necessary capacity of the park of supercondensators is calculated starting from the equation of energy of ( 3)

where E cap Energy requirement of the supercondensator (J), V a Maximum tension of operation of the supercondensator, V b Minimal tension of operation of the supercondensator.

The number of supercondensators series-connected and in parallel can be obtained respectively by ( 4) and ( 5) by using the capacity calculated into (3) on the basis of maximum and minimum value of E cap .

(1)

N PVmin = ∫ 8760 1 P Load ∫ 8760 1 P PVmax (2) N PVmax = max ∫ 24 1 P Load (jour) min ∫ 24 1 P PVmax (jour) (3) C CapReq = 2E cap V 2 a -V 2 b (4) N CapSer ≥ V a V Cap (5) N CapPar ≥ C CapReq N CapSer C Cap

Model of storage of the batteries

The park of batteries is optimized by taking account of PV ield and the capacity of the park of supercapacitors. The proile of charge/discharge of the System with Hybrid Storage is adjusted by using the proile of the supercapacittor to obtain the proile of charge/discharge park of batteries. The number of batteries series-connected is calculated according to (6):

The determination of the number of batteries connected in parallel to increase the capacity of the battery is more complicated. The supplied energy of the batteries and the reduction of capacity due to the current of discharge are taken into account. The energy of the battery accessible in the System to Hybrid Storage depends on the SOC limit deined in the system and is determined by [START_REF] Seeling-Hochmuth | A combined optimisation concept for the design and operation strategy of hybrid-PV energy systems[END_REF]: By this method, the load and the discharge day labourer of the energy of the battery are compared to calculate the number of consecutive days where the battery is insuiciently charged. Consequently, in this situation, the contribution of the photovoltaic ield and the bench of supercondensators to the load is taken into account, which reduces the initial size of the bench of batteries. The initial number of batteries in parallel is calculated from (8) ( 6) International Journal on Interactive Design and Manufacturing (IJIDeM) 1 3

N BattSer = V Sys V Batt (7) E Batt = V Batt Ah cap 1 -SOC limit (8) N BattPar = E max E Batt × N j

Loss of power supply probability (LPSP)

The LPSP is the probability that a loss of supply energy occurs, which means that PV system and the system of storage of energy are not able to feed the load when it is needed. Consequently, it is the fraction of the load which cannot be fed for the period of operation, going from 0 to 1 [START_REF] Semassou | AIDE A LA DECISION POUR LE CHOIX DE SITES ET SYSTEMES ENERGETIQUES ADAPTES AUX BESOINS DU BENIN[END_REF].The LPSP is calculated by (9) when the LPSP is equal to 1, the request is never satisied and a redimensioning of the system is necessary. A LPSP equal to 0 means that the request is satisied constantly, which, according to the criticality of the load, can involve an oversizing of the system with an increased cost.

Economic model based on the concept of cost on cycle of life

The function of optimization consists in minimizing the cost on cycle of life of the system. The lifespan of each component of the system must be given in the calculation of the cost.

According to the chart of the manufacturer, the photovoltaic panels have one lifespan 20 years, which corresponds to the lifespan of the system [9]. The supercondensators are not replaced because they have one lifespan higher than 500,000 cycles [10].To determine the number of replacements of necessary battery throughout the lifespan of the system, a model of lifespan of the battery was built-in in the process of optimization. The model of battery is an energy model of low based on the model used in the program of optimization Homer [11].

The number of replacements necessary for the battery during the lifespan of the system is given by comparing the annual throughput of the battery and the capacity of low of the battery on the cycle of life (10).

Life cycle cost (LCC) includes the value of the initial investment, the cost of replacing the component, the cost of maintenance and repair and the cost of downtime. For a component of the system I, the economic cost of the life cycle (during 25 years) can be expressed by the following equation [START_REF] Khan | Pre-feasibility study of stand-alone hybrid energy systems for applications in Newfoundland[END_REF][START_REF] Navaeefard | Distributed energy resources capacity determination of a hybrid power system in electricity market[END_REF][START_REF] Dehghan | Optimal Sizing of a Hybrid Wind/PV Plant Considering Reliability Indices[END_REF]: with:

(9) LPSP = 1 T T ∫ 1 P PV + P Storage < P Load (10) N RB = min E Bdef E r , CV B LCC I = N i C i + CR i ⋅ K i + CMR i ⋅ PWA ir, R v y i = R v L i -1 , if R v is dividable to L iy i = R v L i , if R v is not dividable to L i
where N i is the number of component i, CI i is the initial investment cost, CR i is the replacement cost, CMR i is the cost of maintenance and repair of component i. PWA and K i are annual and single payment present worth factors, respectively. y i and L i are the numbers of replacements of component i and its lifetime. ir is the real interest rate, R v is the project's lifetime. The total economic cost of the life cycle of the system can then be deduced:

The price of PV panel, the battery and the supercapacitor is 200 euros, 70 euros and 45 euros respectively.

Method of optimization used

A program is conceived to optimize autonomous PV system with hybrid storage with supercapacitor. The function-objective of the program consists in minimizing the economic cost on cycle of life of the system with diferent constraint from volume and probability of loss of power. The program analyzes the possible angles of inclination of PV panel and calculates the angle which makes it possible to receive the maximum of radiation solar during the year. The program of optimization has like step of time the second, able to collect the peaks of power of the load and the luctuations of the proiles of solar radiation.

The program also makes it possible to determine the inancial economy carried out by exploiting the Hybrid System of Storage to supercapacitor and proposes various combinations of panels, batteries and supercapacitors compared to limit LVD ixed with the cost on cycle of life and the corresponding level of cover in energy.

The various models previously benches were implemented under Matlab for the optimal dimensioning of the hybrid system, the calculation of the cost on cycle of life and the illustration of the level of load of the battery for each proile of load. The diagram of optimization developed is as follows (Fig. 5).

The management auditor of energy for the complete system was studied with like imposed constraint the reill of the supercondensator between the points of power. The algorithm of optimization for the control of the simple system of storage and hybrid storage is presented at the following tables (Tables 1,2).

K i = y i ∑ n=1 1 (1 + ir) nl i PWA ir, R v = (1 + ir) R v -1 ir(1 + ir) R v C total = ∑ i LCC i AQ2 AQ3 286

Results and discussions

Peak load without constraints

The results of optimization are presented at Table 3, indicating to the number of photovoltaic panels, of batteries, supercondensators as well as the cost on cycle of life of the system and the consecutive days of insuicient reill undergone by the battery during the year. The optimal ixed angle so that the site produces the maximum annual power of the International Journal on Interactive Design and Manufacturing (IJIDeM)

I Bat =(P Load -P mppt )∕V sys ) I Ucap = 0 I Load = I pv + I Bat + I Ucap (P mppt < P Load )&(Bat SOC > Bat SOCLim )&(Ucap SOC > Ucap SOCLim )&(I storage > Bat lim ) I Bat = Bat lim I Ucap = (P Load -P mppt ∕V sys )-B lim I Load = I pv + I Bat + I Ucap (P mppt < P Load )&(Bat SOC > Bat SOCLim )&(Ucap SOC ≤ Ucap SOCLim )&(I storage ≤ Bat lim ) I Bat = Bat lim I Ucap = (P Load -P mppt ∕V sys )-B lim I Load = I pv + I Bat + I Ucap (P mppt < P Load )&(Bat SOC > Bat SOCLim )&(Ucap SOC ≤ Ucap SOCLim )&(I storage > Bat lim ) I Bat =(P Load -P mppt )∕V sys ) I Ucap = 0 I Load = I pv + I Bat + I Ucap (P mppt < P Load )&(Bat SOC ≤ Bat SOCLim )&(Ucap SOC > Ucap SOCLim ) I Bat = 0 I Ucap =(P Load -P mppt )
1 3
photovoltaic panels of 17 is directed towards the south. No restriction on volumes of components was implemented in the cases described in Table 3. The reduction of the SOC limit of the battery reduces the costs of the system, because more energy is available starting from the bench of batteries.

The optimized results of Table 3 were simulated by using the models of the system. The optimized system returns the battery to 100% of SOC for the majority of the year. There are two periods of weak prolonged solar radiation, one occurring at the beginning of the duration of optimization and the other close to the end. For these prolonged periods, the battery is insuiciently reloaded each day, which involves an additional fall of the load of the battery each day. Figure 6 gives the daily minimal and maximum SOC values of the bench of batteries, for the Hybrid System of Storage over the period of optimization for a SOC limit of 20%.

The return of the bench of batteries to 100% of SOC helps each day to prevent sulphating and the stratiication. Sulphating is caused by the maintenance of the battery in a state of under load for one prolonged period of time. The lead sulphate crystals are usually transformed into lead during the process of load. If the battery is not suiciently reloaded, certain lead sulphate crystals remain on the plates. Finally, soft the lead sulphate crystals become hard and act like an insulating layer preventing the battery from receiving the load. The capacity of the battery is reduced and possibly a cellular short-circuit can occur.

The stratiication caused by an insuicient load of the battery occurs in drowned batteries lead-acid. It occurs when the heavier sulphuric acid in the electrolyte runs at the bottom of the case of the battery by leaving almost pure water on the top. Almost pure water at the top of the case of the battery is not a good conductive of electricity, for this reason, only the lower part of the plates for battery is available to conduct electricity reducing the capacity of the battery.

Oxidation occurs on the higher part of the plates causing corrosion. Corrosion involves permanent damage, reduced the performances of the battery and reduces the lifespan of the battery. Equalization and the gazage can help to prevent the stratiication by agitating the electrolyte [START_REF] Lagorse | Sizing optimization of a stand-Alone street lighting system powered by a hybrid system using fuel Cell, PV and battery[END_REF][START_REF] Cugnet | Peukert's law of a lead acid battery simulated by a mathematical model[END_REF].

The process of optimization was carried out during 1 year by holding account of the lifespan of the components during the calculation of the cost of the system. Consequently, because of the short lifespan of the battery, it proved more economic to increase the number of photovoltaic panels necessary to feed the load for the periods of weak solar radiation, generating an energy of discharge for the periods of high solar radiation. The optimal combination of the bench of batteries and PV system is such as the batteries are reloaded with full load during the major part of the year.

Peak load with constraints of volume

The Hybrid System of Storage with a SOC limit of 20% imposed on the park of batteries was examined for the restrictions of volume on the individual component in the system. The number of PV panels for this case was limited International Journal on Interactive Design and Manufacturing (IJIDeM)

1 3
to 15, that is to say a reduction of 12 panels. The reduction of the number of photovoltaic panels leads to a number of consecutive days of insuicient reill of the battery passing with 8, the number of batteries necessary to satisfy the SOC limit increasing to 40. The reduction of the number of photovoltaic panels increases the need for storage for energy. The battery of the Hybrid System of Storage is increased to allow the food of the load for the long periods of weak solar radiation. When the number of photovoltaic panels is reduced, the quantity of energy produced during the year also decreases.

The increase in the number of batteries necessary and their replacement cost throughout the lifespan of the system make pass the total cost to 17,225 euros. The minimal and maximum daily consumption of the battery during the year can be seen on Fig. 7a.

One second analysis of sensitivity limiting the number of seven batteries in the System of Hybrid Storage with a SOC limit of battery of 20% was also carried out. The limitation of the battery in the Hybrid System of Storage requires the installation of a larger photovoltaic generator, in order to generate the daily requirements in energy for the load, because the bench of batteries is not enough large to provide the load for one prolonged period. In this situation, the number of PV panels is almost doubled to 51, which changes the cost of the system to 12,875 €. The great number of PV panels also increases the quantity of energy produced by the photovoltaic generator during a high solar radiation. The annual minimal daily maximum loading of the battery can be seen on Fig. 7b. The great number of photovoltaic panels brings back the bench of batteries to 100% of SOC during the major part of the year, but the reduction of the number of batteries means that the batteries are discharged on a level from weaker load. The comparison of the cost of the system with life of the two analyses of sensitivity highlights the impact of the battery on the cost on cycle of life because of the number of replacements of necessary batteries for this period.

Peak load with LPSP constraint

The size of the system and the cost on cycle of life can be reduced if the outputs were not critical and if a LPSP higher than zero were taken into account. During the major part of the year, the batteries of the system are maintained on a level of load raised in the system of hybrid storage. If the system were authorized not to feed the load for a short period of time for prolonged periods of weak solar radiation which occur only twice during the year, the park of batteries could be reduced. The adjustment of the LPSP with 0.9132 × 10 -3 in the process of optimization for the system of hybrid storage with a SOC limit of 20% tiny room the cost of the system with 9825 €. This is carried out by reducing the size of the battery. A pressure loss of approximately nine hours was observed for the period of weak solar radiation. The maximum and minimal daily SOC values of the battery for the hybrid system of storage are indicated on Fig. 8. To take a greater limit of LPSP will more reduce the cost on cycle of life of the system to the detriment of a less reliable system.

Conclusions

The results shows that supercapacitors decrease losses in batteries keeping them at a good state of charge thanks to low-level discharges current. In most systems with battery storage, these are sized to support instant powers. In the presence of supercapacitors, this dimensioning will be done according to the average power to drop size of batteries without afecting normal operation of the system. Supercapacitors improve the operation of batteries and so increase their lifespan and then reduce investment costs. Also, the fact that supercapacitors can be completely discharged and stored as is, makes their very interesting use for PV applications.

Results also shows that supercapacitors take in account peaks current and therefore relieve the batteries. In the case where these are placed in outside, the beneical role of supercapacitors is even greater since lead-acid batteries went very badly at low temperatures.

It should be noted that the eiciency of the converters is very low while it is impossible to make a reliable energy balance (losses in energy are greater than the energy provided by supercapacitors). In conclusion, these tests allow to show that supercapacitors allow optimal functioning of batteries, resulting in an increase in their service life.

In conclusion, to achieve the objectives, we explored the optimal design of autonomous photovoltaic solar systems with hybrid storage in accordance with the speciications.

It was necessary to apply knowledge in energy and industrial engineering. It was also necessary to integrate digital simulation into collaborative design environments in order to ind the best candidate for design variables for each piece of equipment in the system. Throughout this approach, the objective of this research is achieved: reduce the risks in the decision-making design, optimize the choice in the preliminary design and develop optimal robust design methods for solar systems with hybrid storage.

International

  Journal on Interactive Design and Manufacturing (IJIDeM)

Fig. 1 a

 1 Fig. 1 a Proile of daily maximum load. b Proile of load of impulse

Fig. 2

 2 Fig. 2 Proile of domestic load

Fig. 4

 4 Fig. 4 Coniguration of autonomous photovoltaic system with hybrid storage

Fig. 5

 5 Fig. 5 Optimization methodology

  (P mppt = P Load )&(Ucap SOC > 100%)|(P mppt = P Load )&(Bat SOC ≤ Bat SOCLim ) I Bat = 0 I Ucap = 0 I Load = I pv + I Bat + I Ucap (P mppt > P Load )&(Ucap SOC < 100%) I Bat = 0 I Ucap =(P Load -P mppt )∕V sys I Load = I pv + I Bat + I Ucap (P mppt > P Load )&(Bat SOC < 100%)&(Ucap SOC ≥ 100%) I Bat =(P Load -P mppt )∕V sys I Ucap = 0 I Load = I pv + I Bat + I Ucap (P mppt > P Load )&(Bat SOC ≥ 100%)&(Ucap SOC ≥ 100%) I Bat = 0 I Ucap = 0 DumpEnergy = P mppt -P Load I Load = I pv + I Bat + I Ucap -DumpEnergy (P mppt < P Load )&(Bat SOC > Bat SOCLim )&(Ucap SOC > Ucap SOCLim )| (I storage ≤ Bat lim )&(Ucap SOC < 100%) I Bat = B lim I Ucap = (P Load -P mppt ∕V sys )-B lim I Load = I pv + I Bat + I Ucap (P mppt < P Load )&(Bat SOC > Bat SOCLim )&(Ucap SOC > Ucap SOCLim )& (I storage ≤ Bat lim )&(Ucap SOC ≥ 100%)

  ∕V sys I Load = I pv + I Bat + I Ucap Else I Bat = 0 I Ucap =-P mppt ∕V sys I Load = I pv + I Bat + I Ucap

Fig. 6

 6 Fig.[START_REF] Borowy | Methodology for optimally sizing the combination of a battery bank and PV array in a wind/ PV hybrid system[END_REF] Maximum/minimum daily SOC of the battery for the hybrid storage system with SOC limit of 20%

Fig. 7

 7 Fig. 7 maximum/minimum daily SOC of the battery (a) with PV ield limited to 15 panels and SOC of 20% (b) with park of batteries limited to 7 batteries and SOC of 20%

Table 1

 1 Controls system of storage of the battery Load = I pv + I Bat (P mppt > P Load )&(Bat SOC ≥ 100%) I Bat = 0 DumpEnergy = P mppt -P Load I Load = I pv + I Bat -DumpEnergy (P mppt > P Load )&(Bat SOC < 100%) I Bat =(P Load -P mppt )∕V sys I Load = I pv + I Bat (P mppt < P Load )&(Bat SOC > Bat SOCLim ) I Bat =(P Load -P mppt )∕V sys I Load = I pv + I Bat Else I Bat = P mppt ∕V sys I Load = 0

	Condition	Action
	(P mppt = P Load )	I Bat = 0
		I

Table 2

 2 Controls hybrid system of storage I

pv = P mppt ∕V sys Condition Action

(P mppt = P Load )&(Ucap SOC < 100%)&(Bat SOC > Bat SOCLim ) I Bat = B lim I Ucap =-B lim I Load = I pv + I Bat + I Ucap

Table 3

 3

	System of hybrid storage versus optimized for	SOC limit	Hybrid storage system			Storage system with only battery
	various limits SOC		N PV	N BAT	Ncap	Cost (€)	Jours	N PV	N BAT	Cost (€)	Days
							InadChar				Inad-
											Char
		20	27	14	1	10,525	7	29	16	11,400	7
		30	29	14	1	10,925	7	33	15	11,850	7
		40	31	15	1	11,675	7	31	19	12,850	7
		50	33	16	1	12,425	7	33	21	13,950	7
		60	33	20	1	13,825	7	33	25	15,350	7
		70	34	26	1	16,125	7	34	32	18,000	7
		80	56	23	1	19,475	4	52	34	22,300	4
		90	83	29	1	26,975	4	75	50	32,500	2


Clément Bonneau

Christophe Josset

Vincent Melot

Bruno Auvity
email: bruno.auvity@univ-nantes.fr

Comprehensive Review of Pure Vapour Condensation Outside of Horizontal Smooth Tubes

Keywords: Condensation, Condenser, Tube bundle, Heat transfer coefficient Condensation, Condenser, Tube bundle, Heat transfer coefficient

published or not. The documents may come

Comprehensive review of pure vapour condensation outside of horizontal smooth tubes

Introduction

General introduction

Thermo-economic optimization of thermodynamic cycles is a topic of current interest [2] due to growing concerns regarding energy use. A key role in these cycles is played by the condensers [3], since they are responsible for most exergy loss [4] - [8]. It is therefore crucial to properly describe this equipment's thermalhydraulic performances. In this context, Naval Group, which is a world leader in naval defence and an innovative player in energy, intends to improve its condenser design tools in order to manufacture more efficient products. Whether they are embedded on military ships, where their size matters, or used in onshore facilities, where their cost matters, it is crucial to properly design the condensers with margins as reasonable as possible.

The most common surface condenser met is the shell-and-tube one, which consists of several circular pipes within a cylindrical shell enclosed by tubes sheets at each end, with the condensation phenomenon occurring on the outside of these tubes (i.e. shellside). Shell-and-tube condensers are found in many applications, since they offer a wide operating range in terms of pressure, fluids and power. Besides, their modularity makes them easy to design, which may be the reason of their success over the 20th century. In spite of that, for economic reasons, their efficiency had to be improved, not only through their designers experience feedback, but through a deeper understanding of condensation phenomena. In the meantime, most condenser designers would use the formulation provided by Nusselt in the early 20th century for condensation of pure stagnant vapour outside a horizontal smooth tube [9], which proved itself to be precise enough to predict condensers performance within an acceptable margin of error. The study of condensation outside of smooth tubes progressively ended in the 1980's to the advantage of 20 enhanced surface tubes. Nowadays, manufacturers are more interested in new technologies such as plate condensers, which offer a wide field of research. However, plate condensers still have limitations that make them unsuitable for specific applications.

For instance, in the electro-nuclear sector, they are nowhere to be found, for they cannot meet the high power requirements. Besides, there are mechanical issues, because there usually is a strong pressure difference between the vapour side and the cooling water side. These mechanical constraints are well handled in shell-and-tube condensers, since it is related to the tube thickness. However, for plate exchangers, the design is a bit more complicated. The size limitation is mainly due to the current manufacturers' equipment. These reasons have led some manufacturers to keep designing shell-and-tube condensers. Therefore, shell-and-tubes condensers are commonly used in refrigeration, air-conditioning 30 and heat pump equipment of medium to large capacity.

As mentioned before, studies of tubes with enhanced surfaces, both inside and outside, have progressively appeared over the 1980s. These new tube geometries have proved themselves quite efficient, but mostly for single-phase flows. As for condensers, the shellside tube enhancement has brought a very small benefit at best [10], since the main thermal resistance comes from the tube-side flow. Fur-35 thermore, when considering industrial aspects, only low finned tube may be considered, because of the several baffle plates electronuclear condensers have: the tube assembly would scrap these fins.

Nowadays, most condenser designs are still based on the standards of the Heat Exchange Institute (HEI) [11], which are based on formulation and data by Orrok (1910, [12]). In the meantime, numerical modelling of industrial condensers appeared in the early 1980s [13] - [25]. These approaches are all 40 based on the resistance summation method, which is often recommended instead of the HEI method [26] [27]. Such a local approach is probably better-suited than a global one, when considering a detailed performance characterisation within an optimisation process. Coupling this method with the porous media approach, which is an homogenisation method used to obtain a sufficient description of the vapour flow, a complete CFD tool is obtained such as the one developed by the present authors [28] 45 [29].

What struck the present authors in above mentioned publications, is the constant change in resistance correlations from an author to another, and even from a publication to another of the same author. Plus, no justification is provided by the authors regarding the reasons that lead to their choice. Since the heat transfer is the core of such modelling, it is of major importance to correctly choose the most appropriate 50 correlation. Indeed, the porous media CFD modelling requires a deep understanding of implemented correlations, otherwise no optimization process can be run.

Therefore the present authors intend to clarify the meaning of the various existing correlations for the heat transfer coefficient in the case of pure vapour condensation outside of horizontal smooth tubes.

Though several authors have already partially reviewed this case [30] [31] [32], they appeared to be 55 incomplete, as they deal with several modes of condensation, and to be sometime inaccurate in the cited literature. The purpose of the present article is therefore to describe how these correlations were obtained and how they are interconnected. It is meant to help CFD engineers, researchers and designers choose the appropriate correlations for their condenser modelling.

The study of condensation in industrial shell-and-tube condensers is complex since it is a combination of elementary phenomena. The methodology followed in the present article is to start from a single tube with numerous assumptions (i.e. Nusselt's work), and to progressively lift them. As shown in Figure 1 where tubes are presented in white and condensate in black, there are two main phenomena that drive the condensation process. On the horizontal axis, there is the increasing impact of the vapour flow 65 around the tube, which causes shear stress on the condensate film, and on the vertical axis, there is the influence of tube bundle effects on local condensation, due to condensate inundation. A few general situations can be highlighted:

• Condensation around a single tube in stagnant vapour (a) ( §3)

• Condensation around a single tube under shear stress, with a laminar condensate film (b) ( §4)

70

• Condensation around a single tube under shear stress, with a turbulent condensate film (c) ( §4)

• Condensation in a bundle of tubes in stagnant vapour, with a laminar condensate film (d) ( §5)

• Condensation in a bundle of tubes in stagnant vapour, with a turbulent condensate film (e) ( §5)

• Condensation in a bundle of tubes under shear stress. This situation lies within the bottom-right corner in Figure 1, where both shear stress and condensate inundation are present at the same time. It is the most complex situation, hence the lack of graphical representation, but also the most common in industrial condensers.

3 Condensation on a single tube with vapour at rest 3.1 Nusselt's theory: The cornerstone A hundred years ago, german professor Wilhelm Nusselt studied analytically the condensation phenomenon on both a vertical plate and a horizontal tube. This pioneering work [9] has been used ever since, mentioned in every publication on surface condensation. As for the horizontal tube study, here are the hypotheses:

1. The surrounding vapour is at rest.

ν l ∂ 2 u ∂r 2 + g sin(θ) = 0 (1)
The pressure gradient term is here neglected since it is much smaller than the body force effect.

95

Then this equation is solved using the following boundary conditions:

     u(r = 0) = 0 (no slip boundary) ∂u ∂r δ
= 0 (no shear stress at vapour-liquid interface)

with δ the local film thickness.

We obtain:

u(r, θ) = g 2ν l sin(θ)r(2δ -r) (2)
and the mean velocity over the thickness of the film is:

u = 1 δ δ 0 u dr = g 3ν l sin(θ)δ 2 (3)
As presented in Figure 2, a mass balance is achieved over a portion of the condensate film, which 100 leads to the following equation:

d(uδρ l ) = j θ • Rdθ ( 4
)
where u is the mean velocity of the condensate over the thickness of the film and j θ is the condensation mass flux in the case of pure conduction defined by:

j θ = λ l ∆T δ∆h v (5)
where ∆T = T sat -T w is the temperature difference between the saturated vapour and the tube wall.

This equation is obtained by a simple energy balance equating the phase change enthalpy (left member)

105
with the conductive heat transfer through the film (right member) :

j θ ∆h v = λ l ∆T δ (6)
Substituting ( 3) and ( 5) into (4), we obtain:

d dθ (δ 3 sin(θ)) = 3ν l Rλ l ∆T ρ l g∆h v • 1 δ (7)
The solution of this differential equation is:

δ = 3ν l Rλ l ∆T ρ l g∆h v • 4 3(sin(θ)) 4/3 θ 0 (sin(ω)) 1/3 dω 1/4 (8)
Under the assumption that heat is transferred by pure conduction, the local heat transfer coefficient h is:

110 h = λ δ (9)
Finally, integrating h over the tube, we obtain:

h = 1 π π 0 hdθ (10)
= 0.728

ρ l gλ 3 l ∆h v ν l D o ∆T 1/4 (11)
This heat transfer coefficient is often written with the constant equal to 0.725. This is due to a lack of precision in Nusselt's calculations, who probably used a handmade Riemann integral back in the day.

The Nusselt number can therefore be expressed as:

N u = 0.728 ρ 2 l gD 3 o ∆h v µ l λ l ∆T 1/4 (12)
This work is the cornerstone of all scientific production on filmwise condensation outside horizon-115 tal tubes. The advantage of this method is that it is purely analytical, though limited by certain assumptions. However, it is not limited to a specific fluid, which is often the case with experimental correlations.

Pressure gradient effect

Another form of this equation is commonly found, namely :

120 N u = 0.728 ρ l (ρ l -ρ v )gD 3 o ∆h v µ l λ l ∆T 1/4 (13)
This comes from the mechanical equilibrium equation which takes into account the pressure gradient, which is hydrostatic. Therefore the equation ( 1) is slightly changed :

µ l ∂ 2 u ∂r 2 + (ρ l -ρ v )g • sin(θ) = 0 (14)
It appears that the Nusselt number obtained from ( 13) is equal to one obtained from (12) when

ρ l ρ v .
Belghazi et al. [33] confronted their experimental results for condensation of HFC134a on tubes 125 of the upper row of a bundle with equation ( 13) and obtained discrepancies of about 10%. So did Fernández-Seara et al. [10] who obtained similar discrepancies with ammonia.

Enthalpy convection effect

In the previous analysis, it was assumed that heat was only transmitted by conduction through the condensate film. In the heat balance, the impact of the thermal capacity c p is neglected. In 1952,

Bromley [34] included this thermal capacity in Nusselt's work, since he assumed that heat couldn't only be transmitted by conduction, which case provided with a linear temperature profile (red plain line

Figure 3), as if the condensate was still. According to Bromley, for a given thickness of the condensate film and temperature difference, the temperature profile should be curved (red dotted line Figure 3).

In the case of negligible viscous dissipation and negligible compressibility effects (Bejan, 2004 [35]), 135 10 the stationary energy equation within the condensate film is:

ρ l c p ∇ • T V = λ l ∇ • ∇T (15)
which becomes

∇ • ρ l c p T V -λ l ∇T = 0 (16)
Then, using the divergence theorem over the domain shown in Figure 4 and neglecting the orthoradial component of the temperature gradient:

-λ l ∂T ∂r r=0

Rdθ + λ l ∂T ∂r r=δ

Rdθ

+ H f -(H f + dH f ) = 0 ( 17
)
with H f the enthalpy flow rate per unit length entering the shaded domain:

140 H f = δ 0 ρ l uc p (T -T sat ) dr (18)
The incoming heat flux on the right side of the domain comes from the condensation mass rate:

λ l ∂T ∂r r=δ = j θ • ∆h v (19)
In order to obtain the proper temperature profile, Bromley solved the same equations as Nusselt did, except that he corrected (6) by using the above equations. This calculation requires the temperature gradient at the wall, which is not known a priori. Therefore, Bromley proceeded with an iterative procedure and initialised with a linear temperature profile. At the end of the first iteration, he obtained:

145 N u = 0.728 1 + 3 8 Ja 1 + 7 30 Ja 3/4 ρ l (ρ l -ρ v )gD 3 o ∆h v µ l λ l ∆T 1/4 ( 20
)
Where Ja is the Jakob number. For c p = 0, (20) is equal to (13). Then, Bromley repeated the same procedure, starting from the new temperature gradient at the wall. At the end of the second iteration, he obtained the same equation as above with 1 + 3 8 Ja 1 + 7 30 Ja 3/4 replaced by:

1 + 1 2 - 1 8
1 + 0.052322 • Ja 1 + 0.233333 • Ja • Ja 1 + 1 3 - 1 10
1 + 0.051786 • Ja 1 + 0.233333 • Ja • Ja 3/4 (21)
This second approximation is really close to the first one, therefore no further calculation was achieved. Bromley proposed a simpler correction:

150 √ 1 + 0.4 • Ja ( 22
)
This correction is often applied to the phase change enthalpy, namely : 1 for a wide range of Ja numbers. The simplified version is a good approximation for Ja below unity, but differs when above. However, for steam condensation, Ja rarely exceeds 0.1.

∆h v = ∆h v (1 + 0.4 • Ja) 2 (23)
In 1956, Rohsenow [36] did the same analysis than Bromley [34] for a vertical plate but took into 155 account the effect of the crossflow within the condensate film. The force balance for the vertical plate is similar to the one from the tube and therefore the velocity is:

w(x, z) = g(ρ l -ρ v ) µ l xδ(z) - x 2 2 (24)
For the red control volume in Figure 5, the mass flowrate d ṁ across the face at x is:

d ṁ(x) = x 0 ρ l w(x, z + dz)dx - x 0 ρ l w(x, z)dx = gρ l (ρ l -ρ v ) µ l x 2 2 δ ( 25
)
Then the heat balance at a distance x from the wall:

λ l ∂T ∂x = j θ ∆h v + c p ρ l w(x)T (x) + c p d ṁ(x)T (x) (26)
Repeating the same procedure as Bromley's, the obtained simplified correction after a few iterations is:

∆h v = ∆h v (1 + 0.68 • Ja) (27)
It is close to Bromley's correction, which is more obvious when ( 23) is expanded into Taylor series for small values of Ja:

∆h v = ∆h v (1 + 0.8 • Ja) (28)
There is no apparent restriction in applying Rohsenow's correction for a horizontal tube.

In 1959, a boundary-layer analysis for laminar film condensation was performed by Sparrow & Gregg 165 for both a vertical plate [37] and a horizontal tube [38]. Their analysis took into account both the inertia forces and energy convection. They considered the following system of equations:

∂u ∂x + ∂v ∂y = 0 (Mass) (29)
u ∂u ∂x + v ∂u ∂y = g • sin x r + ν l ∂ 2 u ∂y 2 (Momentum) (30)
u ∂T ∂x + v ∂T ∂y = λ l ρ l c p ∂ 2 T ∂y 2 (Energy) (31)
The system of coordinates is described in Figure 6.

Then they reduced these partial differential equations to ordinary differential equations, which have from 0 to 1. However, no formulation of the Nusselt number depending on P r and Ja could be expressed.

For small value of Ja, the Nusselt number tend to:

N u = 0.733 ρ 2 l gD 3 o ∆h v µ l λ l ∆T 1/4 (32)
which is really close to (12). Sadasivan & Lienhard (1987 [39]) noticed that Bromley and Rohsenow's developments involved the implicit assumption that the P r was infinite. However, Sparrow & Gregg demonstrated that the Nusselt 175 number depended upon P r and Ja, and that for high values of P r, their solution was close to Nusselt's solution with Rohsenow's correction. So the authors assumed that the general solution was of the form:

∆h v = ∆h v (1 + C(P r) • Ja) (33)
where C was a parameter to determine.

For P r ranging from 0.6 to 1000 and Ja below 0.8, Sadasivan & Lienhard obtained the following formulation using a best curve fitting method:

180 ∆h v = ∆h v 1 + 0.683 - 0.228 P r l Ja (34)

Surface temperature variation

All previous work assumed an isothermal surface at the tube wall for practical reasons. Nevertheless, due to an uneven film thickness around the tube, the wall temperature is higher at the top of the tube, where the film is thinner, and lower at the bottom, where the film thickness becomes near infinite. This phenomenon has been studied by Memory & Rose [40], who assumed a cosine temperature distribution:

185 T w = A • cos(θ)∆T + T w (35)
Where

A ∈ [0; 1].
This assumption is based on the conclusion of Lee et al. and obtained the Nusselt constant equal to 0.7280 quite surprisingly. These first four significant digits remained constant for all values of A.

(3) into (37):

d dθ (δ 3 sin(θ)) = 3µ l D o q 2ρ 2 l g∆h v (38)
The solution of this equation is:

δ = 3µ l D o q 2ρ 2 l g∆h v • 1 sinc(θ) 1/3 (39)
Then the mean heat transfer coefficient is obtained: (40) and the corresponding Nusselt number:

h = 0.693 λ 3 l ρ 2 l g∆h v µ l D o q 1/3
N u = 0.693 ρ 2 l gD 2 o ∆h v µ l q 1/3 (41)
The analysis of Butterworth [43] regarding this work is quite confusing, since it suggests that the 200 only difference between ( 12) and ( 41) is the coefficient changing from 0.728 to 0.693, whereas the whole expression is different.

As mentioned by the authors, the assumption of uniform wall temperature makes sense when the cooling side is the most resistive, but it is no longer appropriate when both sides are of closer magnitudes.

4 Condensation on a single tube under vapour flow 205 In the previous section, the vapour surrounding the tube was at rest and the condensate film was set in motion by gravity. The following section will present the phenomenology of filmwise condensation under vapour crossflow. Indeed, for high vapour velocities, the condensate film undergoes a shear stress at liquid-vapour interface. This shear stress enhances the heat transfer by locally thinning the liquid film. For a vertical 210 downflow as represented in Figure 7, the film is thinner on the upper half while it is thicker on the lower half. The thickening is caused by the separation of the vapour boundary-layer after the separation point, which is located at an angle usually comprised between 80 • and 180 • from the top of the tube.

At this very point, the shear stress at the liquid-vapour interface changes sign, hence the recirculation flow as presented in Figure 8. Below this separation point, the vapour flows in the direction opposite to 215 the gravity, which causes the thickening and considerably lowers the heat transfer. This separation is caused by the vapour pressure gradient, which is more important when the vapour is flowing, due to the dynamic component of the pressure. This gradient is positive on the forward half, which tend to thin the film, while it is negative on the rear half, provoking the opposite effect.

Approach of Sugawara et al. (1956) 220

Though condensation of a vapour flow on a single horizontal tube has been treated by Fuks [44] in the USSR, the first semi-analytical study was conducted by Sugawara et al. [45] from Japan. They kept most of Nusselt assumptions, except the fact that the vapour is not at rest. Therefore, they obtained the same equations ( 1) & ( 4), but the momentum equation is solved using different boundary conditions, namely:

225      u(r = 0) = 0 (no slip boundary) µ l ∂u ∂r δ = 1 2 f ρ v U 2 ∞
(shear stress at vapour-liquid interface) with f the dimensionless friction factor.

Consequently:

u(r, θ) = g 2ν l sin(θ)r(2δ -r) + f ρ v U 2 ∞ 2µ l r (42)
and the mean velocity over the thickness of the film is:

u = g 3ν l sin(θ)δ 2 + f ρ v U 2 ∞ 4µ l δ (43)
Substituting ( 43) into (4):

ρ l g 3 d dθ δ 3 sin(θ) + ρ v U 2 ∞ 4 d dθ f δ 2 = ν l Rλ l ∆T ∆h v • 1 δ (44)
The friction factor is estimated by using experimental results for a single phase flow around a cylinder.

over the rear half of the tube. Therefore, the results may be too optimistic.

Approach of Shekriladze & Gomelauri (1966)

In 1966, Shekriladze & Gomelauri [46] noticed that previous work ([44] [45] [47]) used the assumption that the shear stress at the liquid-vapour interface was the same as at a dry tube surface without condensation, thus ignoring the momentum transfer caused by the mass of the condensing vapour.

The authors took it into account in their study, under similar assumptions to Nusselt's. However, they assumed that the tube was in a vertical downward potential flow of vapour and inertia forces were neglected. In a potential flow with a velocity far from the tube U ∞ , the velocity varies from 0 at the front and back stagnation points to 2U ∞ on the sides:

U (θ) = 2U ∞ sin(θ) (46)
The following set of equations describe the condensate film flow:

µ l ∂ 2 u ∂r 2 = 0 (Momentum) (47) d(uδρ l ) = j θ • Rdθ (Heat balance) (48)
with the boundary conditions:

     u(r = 0) = 0 (no slip boundary) µ l ∂u ∂r δ
= τ δ (shear stress at vapour-liquid interface)

with the local shear stress defined by Shekriladze & Gomelauri as:

τ δ = j θ (U (θ) -U lv ) (49)
with U lv U (θ) the velocity at liquid-vapour interface (neglected).

The local heat transfer coefficient obtained from this set of equations is:

h(θ) = sin(θ) 1 -cos(θ) λ 2 l ρ l U ∞ µ l D o (50)
Hence the mean heat transfer coefficient:

250 h = 2 √ 2 π λ 2 l ρ l U ∞ µ l D o (51)
Introducing the two-phase flow Reynolds number:

Re = ρ l DU ∞ µ l (52)
The corresponding Nusselt number is:

N u SH = 0.900 Re 1/2 (53)
Being unable to analytically solve (47) with added gravitational term, the authors used an asymptotic model, which formulation approaches their previous result for high vapour velocities (i.e. N u SH ) and the basic Nusselt formulation for low or nil vapour velocity (i.e. N u GR ), namely:

255 N u = N u 2 SH 2 + N u 4 SH 4 + N u 4 GR 1/2 1/2 (54)
Introducing the dimensionless number F :

F = P r F r • Ja (55)
The previous equation may be expressed as a function of F :

N u = 0.637 1 + (1 + 1.68F ) 1/2 1/2 Re 1/2 (56)
This dimensionless number describes if the condensate film is subject to either vapour drag (F -→ 0)

or gravity (F -→ ∞).

In (56), the original coefficient 0.725 from Nusselt was used, instead of 0.728.

260

However, they did not take into account the effect of the pressure gradient, which results in a too optimistic Nusselt number. Well aware of this phenomenon, Shekriladze & Gomelauri proposed a correction of (56). According to the authors, the boundary-layer separation point only happens beyond 82 • , or 65% of the heat transfer takes place on the tube surface comprised between 0 • and 82 • .

Therefore, they decided to neglect the heat transfer on the surface lying beyond this angle, where the boundary-layer may be separated, thus obtaining a 35% lower Nusselt number:

N u = 0.414 1 + (1 + 1.68F ) 1/2 1/2 Re 1/2 (57)
However, Butterworth [48] noticed that reducing the whole Nusselt number by 35% also affect the solution for stagnant vapour case and thus proposed the following correction:

N u = (0.65N u SH ) 2 2 + (0.65N u SH ) 4 4 + N u 4 GR 1/2 1/2 (58)
which results in:

N u = 0.414 1 + (1 + 9.42F ) 1/2 1/2 Re 1/2 (59)
This is the conservative version of (56), which safely underestimates the heat transfer coefficient. for equations ( 12), ( 56) and ( 59) as a function of F .

For F -→ ∞, all formulations approach the Nusselt solution.

It should be noted that the above equations are particularly conservative, given the fact that behind the separation point the heat transfer is assumed to be nil while it is just lowered. Furthermore, the which is less than half of the tube surface. Depending on the actual angle of separation, the Nusselt number is likely to lie between equations ( 56) and (59).

Approach of Fujii et al. (1972)

Unsatisfied with previous work, Japanese searchers Fujii et al. [1] extended the analysis to the vapour boundary-layer. Therefore, no assumption is made regarding the shear stress at the interface. Beyond 280 the boundary layer, the vapour flow is assumed to be a vertical downward potential flow. They obtained the following set of equations for the condensate film:

1 r ∂u ∂θ + ∂v ∂r = 0 (Continuity) ( 60
)
ν l ∂ 2 u ∂r 2 + g • sin(θ) = 0 (Momentum) (61)
d(uδρ l ) = j θ • Rdθ (Heat balance) (62)
and for the vapour boundary-layer:

1 r ∂U ∂θ + ∂V ∂r = 0 (Continuity) (63)
U r

∂U ∂θ + V ∂V ∂r = ν v ∂ 2 U ∂r 2 + 2U 2 ∞ r sin(2θ) (Momentum) (64)
with the boundary conditions:

             u(r = 0) = 0 (no slip boundary) µ l ∂u ∂r δ = µ v ∂U ∂r δ
(shear stress at vapour-liquid interface)

ρ l u R ∂δ ∂θ -v δ = ρ v U R ∂δ ∂θ -V δ
(mass transfer at vapour-liquid interface)

For convenience of solving, the authors decided to neglect the pressure term in (61). They solved

N u SH = 0.90 (RH) -1/3 Re 1/2 (65)
A general solution for all RH was obtained using an asymptotic model of the form:

(N u 3 1 + N u 3 2 ) 1/3 (66)
Introducing the coefficient χ:

χ = 0.90 1 + 1 RH 1/3 (67)
The obtained Nusselt number is:

295 N u SH = χ Re 1/2 (68)
The results for small oncoming vapour velocity tended to Nusselt's equation (12) with the original coefficient of 0.725. Then to connect N u SH (68) and N u GR (12), another asymptotic model was used, namely:

N u = N u 4 SH + N u 4 GR 1/4 (69)
which expressed as a function of F :

N u = χ 4 + 0.276F 1/4 Re 1/2 (70)
Fujii et al. confronted this correlation to experimental results of condensation of steam on a single 300 horizontal tube [1] and on banks of horizontal tubes [49]. They noticed a fair agreement, except for the case of the tube bank with in-line arrangement, where heat transfer was about 20% lower than expected.

This discrepancy may be justified by the specific flow pattern. Therefore, the authors added a coefficient that equals either 1 or 0.8.

However, Fujii and coworkers made the assumption that the vapour flow outside the boundary-layer was a potential one, which causes ∂U ∂r δ to always be positive, therefore preventing the boundary-layer from separating.

Their solution was conservatively modified by Lee & Rose [50] from UK, whose analysis neglects heat transfer behind the calculated separation point. Though little information is available on their method, the same Lee mentioned in his PhD thesis [51] the following equation, which is due to Prandtl:

310 -v δ U ∞ Re v = 4.36 -cos(θ c ) (71)
where θ c > π/2 is the critical angle at which separation occurs.

Then, using this angle, they correlated again the value of N u SH and obtained the following coefficient χ:

χ = 0.88 1 + 0.74 RH 1/3 (72)
Concerning N u GR , the Nusselt equation was used with the more precise coefficient 0.728, hence the mean Nusselt number:

315 N u = χ 4 + 0.281F 1/4 Re 1/2 (73)
Equations ( 70) and ( 73

Approach of Rose (1984)

In all previous studies, the pressure gradient within the condensate film, arising from the vapour flow 320 around the tube, has been omitted in the momentum equation. The effect of this pressure gradient upon the condensate film has been studied by Rose [52], whose analysis was based upon the same assumptions as Shekriladze & Gomelauri [46], except he took into account the body force and pressure gradient within the condensate film. Though the asymptotic value of the surface shear stress was known to be inaccurate, it was adopted owing to its simplicity and the fact that Rose's aim was to investigate 325 the effect of the pressure term. The following set of equations was solved:

µ l ∂ 2 u ∂r 2 + ρ l g • sin(θ) - 1 R dp dθ = 0 (Momentum) (74) d(uδρ l ) = j θ • Rdθ (Heat balance) (75) 22
with the boundary-layer conditions:

     u(r = 0) = 0 (no slip boundary) µ l ∂u ∂r δ
= τ δ (shear stress at vapour-liquid interface) Equation ( 74) may be integrated to give:

u(r, θ) = 1 µ l τ δ r -ρ l g • sin(θ) - 1 R dp dθ r 2 2 -δr (76)
Assuming a potential flow outside the vapour boundary-layer yields:

τ δ = 2j θ U ∞ sin(θ) ( 77
)
dp dθ = -2ρ v U 2 ∞ sin(2θ) (78)
Injecting ( 77) and ( 78) into (76) gives:

u(r, θ) = 1 µ l 2j θ U ∞ sin(θ)r -ρ l g • sin(θ) + 2ρ v U 2 ∞ R sin(2θ) r 2 2 -δr (79)
Then the mean velocity over the thickness of the film:

u = ρ l g • sin(θ) + 2ρ v U 2 ∞ R sin(2θ) δ 2 3µ l + j θ U ∞ δ µ l sin(θ) (80)
Injecting the expression of u and the definition of j θ in (75):

1 δ λ l ∆T ∆h v = ρ l µ l R d dθ ρ l g • sin(θ) + 2ρ v U 2 ∞ R sin(2θ) δ 3 3 + λ l ∆T U ∞ ∆h v δ sin(θ) (81)
Two dimensionless groups can be formed: the dimensionless film thickness δ *

δ * = δ U ∞ Rν l 1/2 (82)
and the dimensionless parameter P which relates to the inclusion of the pressure gradient

335 P = ρ v ∆h v ν l λ l ∆T (83) = ρ v ρ l H -1 (84) = ρ v ρ l • P r Ja (85)
Re-arranging (81):

1 δ * = d dθ F 2 • sin(θ) + 2P • sin(2θ) δ * 3 3 + δ * sin(θ) (86)
with the boundary condition due to the symmetry at the top of the tube:

dδ * dθ = 0 at θ = 0 (87)
The last term inside the brackets in (86) results from the shear stress at vapour-liquid interface, while the second one relates to the inclusion of the pressure gradient in the analysis. If both of these terms are omitted, then equation (86) reduces to the simple equation ( 7) from Nusselt [9].

340

Then, Rose considered the possibility that the velocity gradient at the wall may be nil or negative, which would lead to a separation of the condensate boundary-layer. Therefore the surface shear stress upstream this separation point is negative:

∂u ∂r r=0 ≤ 0 (88)
which becomes:

1 + δ * 2 2 F 2 + 4P • cos(θ c ) ≤ 0 ( 89
)
This inequality is satisfied for θ > π/2 and P > F/8, which means that a solution of (86) only 345 exists for θ ∈ [0; θ c ], unless θ c = π, in which case the equation may be solved over the whole tube.

Then the authors numerically solved (86) for F = 10 3 , 10 2 , 10, 1, 10 -1 , 10 -2 , 10 -3 , 0 and for P = 0, 0.001, 0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1, 2, 5, 10.

The present authors have solved this differential equation (86) using Runge-Kutta fourth-order method.

P = 0

In this case, the pressure gradient is omitted, which corresponds to the analysis of Shekriladze & Gomelauri [46]. Rose correlated his results and obtained:

N u = 0.9(1 + (RH) -1 ) 1/3 + 0.728F 1/2 (1 + 3.44F 1/2 + F ) 1/4 Re 1/2 (90) 4.4.2 0 < P ≤ F/8
When the pressure gradient is taken into account (i.e. P > 0), Figure 12 shows that the condensate film is thinner on the upper part (θ < π/2) while it is thicker on the lower half. Rose noted only small differences for the mean Nusselt number for P = 0 and for P > 0. Even for the limiting cases were P was close to F/8, the maximum discrepancy was around 5%. Therefore, Rose advocated the use of equation (90).

P > F/8

In this case, no solution can be obtained beyond the separation point (i.e. beyond θ c ). As shown in Figure 13 and Figure 14, the closer is θ c from π/2, the thinner is the condensate film on the upper part.

The author's suggestion was to obtained an accurate solution over the upper part of the tube, where no boundary-layer separation can occur, and to neglect heat transfer over the lower half, thus obtaining a conservative solution.

365

Let N u π/2 be the mean Nusselt number for the upper half of the tube. On one hand, when the pressure gradient is omitted (i.e. P = 0), N u π/2 is given by:

N u π/2,P =0 = 1.273 + 0.866F 1/2 (1 + 3.51F 0.53 + F ) 1/4 Re 1/2 (91)
On the other hand, when it is included, for the case of large oncoming vapour velocity (i.e. F -→ 0), N u π/2 is given by:

N u π/2,P =0 = 1.273(1 + 1.81P ) 0.209 (1 + (RH) -1 ) 1/3 Re 1/2 (92)
Combining equations ( 91) and ( 92), for the upper half of the tube:

370 N u π/2 = 1.273(1 + 1.81P ) 0.209 (1 + (RH) -1 ) 1/3 + 0.866F 1/2 (1 + 3.51F 0.53 + F ) 1/4 Re 1/2 (93)
Finally, equation ( 93) is then applied to whole surface of the tube, without considering any other heat transfer. Thus the mean Nusselt number should be divided by 2, but in order to maintain the convergence towards Nusselt's solution for F -→ ∞, the coefficient before F at the numerator is taken equal to 0.728:

N u = 0.636(1 + 1.81P ) 0.209 (1 + (RH) -1 ) 1/3 + 0.728F 1/2 (1 + 3.51F 0.53 + F ) 1/4 Re 1/2 (94)
Rose concluded his article by stating that the pressure gradient effects should be more notable for low 375 vapour flow for refrigerants, due to their higher vapour density, than for steam at equivalent operating conditions. However, the pressure gradient would have a significant effect on steam condensation at high pressures.

Figure (15) shows the pressure gradient effect on N u/ Re

1/2
for different values of P . The noticeable discontinuities in Rose curves (dotted lines) are located at F = 8P , and are due to the switch from 380 equation ( 90) to (94).

In the same paper, Rose confronted his analytical results against experimental ones obtained from the literature for steam, R113 and R21. The results were in really good agreement, which supports the use these formulations.

Approach of Homescu & Panday (1999)

385

In previous cited literature, both vapour and liquid phases were solved assuming laminar flows. In 1999,

Homescu & Panday [53] studied the influence of turbulence in the case of forced convection condensation on a horizontal tube. They also retained in their analysis the pressure gradient, inertia and enthalpy convection terms. Lacking information on local flow structure, the authors assumed that the flow is turbulent all around the tube. Using the above mentioned assumptions, they obtained the following set 390 of equations for the condensate film, in the coordinate system defined in Figure 6:

∂u ∂x + ∂v ∂y = 0 (Continuity) (95) ρ l u ∂u ∂x + v ∂u ∂y = g • sin θ - ∂p ∂x + ∂ ∂y (µ + µ t ) l ∂u ∂y (Momentum) (96) ρ l c p u ∂T ∂x + v ∂T ∂y = ∂ ∂y (k + k t ) l ∂T ∂y (Energy) (97)
and for the vapour boundary-layer:

∂U ∂x + ∂V ∂y = 0 (Continuity) (98) ρ v U ∂U ∂x + V ∂U ∂y = - ∂p ∂x + ∂ ∂y (µ + µ t ) v ∂U ∂y (Momentum) (99) (100)
with the boundary conditions:

                                   u(y = 0) = 0 (no slip boundary) T (y = 0) = T w (isothermal wall) u(y = δ) = U (y = δ) (no slip interface) T (y = δ) = T sat (thermal continuity) (µ + µ t ) l ∂u ∂y δ = (µ + µ t ) v ∂U ∂y δ
(shear stress at vapour-liquid interface)

ρ l u ∂δ ∂x -v δ = ρ v U ∂δ ∂x -V δ
(momentum transfer at vapour-liquid interface)

Regarding the turbulence modelling, the authors used the mixing length concept, where µ t is defined by:

µ t = ρL m ∂u ∂y (101)
where L m is the mixing length.

The authors tested several combinations of models and concluded that the combination of Pletcher's model for the vapour phase and Kato's model best fitted the experimental results. These models are well enough described in the article [53].

Then, they numerically solved this set of equations using a finite difference scheme. To represent 400 their numerical results, the following equation was obtained by the authors:

N u t = 0.291 0.75 1 + (RH) -1 1/3 + 0.25F 1/4 + (1 + 0.8F ) 1/2 (0.25F 1/2 + 1.75F ) 1/4 3/2 Re 1/2 (102)
This Nusselt number describes both the cases of laminar and turbulent condensation. Its evolution is shown in Figure 16 along with Fujii et al. formulation. Unlike any other correlation, this one does not approach Nusselt's curve for F -→ ∞. It approaches the curve 0.173F 3/8 instead of 0.728F

Synthesis

This section sheds light on the most prominent publications related to condensation under shear stress and provides a better understanding of physical phenomena taken into account in each work. The vapour flow tends to thin the condensate film and therefore decrease the thermal resistance of this film. From this review, it would seem appropriate to use the formulation of Rose [52] (equations ( 90) and ( 94)) for it is the most complete approach, considering the phenomena taken into account, and the fact that it remains valid for stagnant vapour. However, it should be noted that it results from an interpolation, and therefore should not be used if parameters are out of the initially considered domain (F < 10 3 and P < 10). Nevertheless, this domain may be expanded by solving the differential equation for the range of interest.

As the vapour flows deeper inside a tube bundle, its velocity decreases, which decreases the condensation rate. Therefore, the tubes on the outer layer of the bundle have the best condensation rate,

435

which is something a condenser designer will benefit by modifying the tubes layout. This is known as geometry effects.

Condensation on a bundle of tubes

As soon as a condenser contains several tubes in a vertical bank, the inundation phenomenon must be accounted for. It is encountered in every single industrial condenser with horizontal tubes. This 440 phenomenon is quite simple, since it only results from condensate streaming down from the bottom of a tube to the top of the tube beneath. The condensate is mostly set in motion by gravity and therefore keeps flooding the tubes up until it reaches the condenser well at the bottom of the device. Then the surrounding condensate film of flooded tubes is thicker and the thermal resistance is increased. This also leads to a higher subcooling of the condensate, which also decreases the global performance of the device. Thus inundation is an undesired phenomenon, that must be dealt with when designing a condenser.

Moreover, this phenomenon depends on many factors, such as the tube bundle geometry, the tube layout and the spacing. The more tubes a vertical bank contains, the more flooded are the lower tubes (see Figure 18(a)). Besides, if the tubes are staggered, then the condensate may not fall directly on the 450 tube beneath as shown in Figure 18(b), but on an intermediate tube, which is named lateral drainage.

When the flow rate is high enough, the drainage splashes on the top of the tube (see Figure 18(c)), and therefore modifies the thickness of the film, which becomes "wavy". This thickness distribution actually gives a better heat transfer coefficient than an even distribution, which partially compensates the thickening of the film.

Correlations based on row number

The first analysis of inundation found in literature dates back from 1949 with Jakob [56]. Starting from Nusselt's theory [9], he adapted the analysis for a vertical bank of isothermal tubes. Let γ n be the local condensation mass rate (on the n-th tube) and Γ n be the global condensation mass rate (over the n tubes) defined as:

460 γ n = πD o h n ∆T ∆h v (105)
Γ n = nπD o h n ∆T ∆h v (106)
Equation ( 12) may be written as a function of γ 1 , which is the condensation mass rate of the top tube (or first tube row):

h 1 λ l µ 2 l ρ 2 l g 1/3 = 1.523 4γ 1 µ l -1/3 (107)
This expression is sometimes used to define a relation between a Nusselt number (left side of the equation) and the condensate Reynolds number (right side of the equation) [43].

465

Assuming that the condensate film falls vertically from one tube to another as a continuous sheet (see Figure 19(a)), then the mean heat transfer coefficient over n tubes is:

h n λ l µ 2 l ρ 2 l g 1/3 = 1.523 4Γ n µ l -1/3 (108)
Then, the inundation factor may be obtained. It represents the decrease of heat transfer due to inundation:

h n h 1 = n -1/4 (109)
and for the n-th tube:

h n h 1 = n 3/4 -(n -1) 3/4 (110)
Thus, under Jakob assumptions, the heat transfer coefficient decreases with the fourth root of 1/n.

For example, for a vertical bank of 16 tubes, the mean heat transfer coefficient is only 50% of the heat transfer coefficient of the top tube. This formulation is sometimes associated to Nusselt himself, though nothing seems to prove that he made these calculations.

According to Kern [57], the condensate is more likely to fall as droplets or columns, rather than as 475 sheets. Such a drainage improves the splashing, and therefore the heat transfer coefficient. The author suggested to replace n by n 2/3 in equation ( 109):

h n h 1 = n -1/6 (111)
and for the n-th tube:

h n h 1 = n 5/6 -(n -1) 5/6 ( 112
)
This formulation is often cited for the design of condensers for it gives pertinent results, while Jakob formulation is too conservative. These two formulations are meant for vertical banks of tubes, where tubes are numbered as shown in Figures 20(a However, in numerous tube bundles with a triangular layout, lateral drainage may be observed depending on the pitch to diameter ratio p/D. This phenomenon has been studied by Eissenberg & Noritake [58] in 1970. When pure lateral drainage is occurring, as shown in Figure 21, the authors obtained the following theoretical result:

h n h 1 lat = 0.6 + 0.4353 n -1/5 (113)
Though no explanation could be found about how this expression was obtained, the present authors obtained similar results with the following calculations. Let N u [0;π/2] and N u [π/2;π] be respectfully the mean Nusselt numbers over the upper half and the lower half obtained from (10) with different intervals of integration:

490 N u [0;π/2] = 0.866 F 1/4 (114)
N u [π/2;π] = 0.590 F 1/4 (115)

Assuming that the condensate drains over the following tube at θ = π/2, then the Nusselt number for the upper half of this second tube remains unchanged, while the Nusselt number for the lower half decreases due to inundation, thus:

N u n N u 1 = 1 2 × 0.866 F 1/4 0.728 F 1/4 + 1 2 × 0.590 F 1/4 n -1/5 0.728 F 1/4 (116)
= 0.5949 + 0.4052 n -1/5 (117)

which is close to (113), especially when considering that Eissenberg & Noritake used the original con-495 stants in their calculations instead of the more precise ones. The coefficient 1/5 for the inundation will be discussed later.

However, according to Figure 20(c), n should actually be n/2 since the condensate drains over an intermediate tube:

h n h 1 lat = 0.5949 + 0.4655 n -1/5 (118)
which is still close to (113).

500

Finally, it was assumed that both left and right halves of the tube were flooded from π/2 to π. It would seem more appropriate to consider that the pattern shown in Figure 21 repeats itself for every column. In this case, only a quarter of the tube is flooded, hence the following inundation factor:

h n h 1 lat = 0.7974 + 0.2327 n -1/5 (119)
For these formulations, it should noted that the inundation factors is not equal to 1 for n = 1 and should therefore be set to this value for the first row.

505

Then Eissenberg & Noritake suggested that a classical vertical drainage would give the following inundation factor:

h n h 1 vert = n -1/5 ( 120
)
This expression is more conservative than Kern's, but less than Jakob's. Though they claimed it was obtained experimentally, no data could be found. However, the inundation factor used in (113) probably originated from here.

510

Then the authors defined a spacing parameter F d , which accounts for the tube layout and ratio p t /D o . Though brief, this parameter is defined in Table 3. Then, this parameter will take into account the proportion of lateral drainage and vertical drainage, using a convex combination: Figure 22 shows the local inundation factors for the above mentioned correlations over 40 tubes. For every correlation, a sharp decrease is noticeable over the first five tubes.

h n h 1 = F d h n h 1 lat + (1 -F d ) h n h 1 vert (121) = 0.6F d + (1 -0.5647F d ) • n -1/
Finally, Short & Brown [59] have obtained experimental results with steam and Freon-11 on a vertical 520 bank of 20 tubes. Their results tend to prove that inundation is negligible due to the counteracting mixing action of the drainage. They concluded that "The average condensate film heat transfer coefficient for a bank of twenty tubes is, to a good approximation, very nearly equal to the value predicted by the Nusselt theory for the top tube in the bank". Butterworth [60] [43] had a different interpretation of these results and claimed that the following inundation factor could be deduced:

h n h 1 = 1.24n -1/4 (124)
and for the n-th tube:

h n h 1 = 1.24 n 4/5 -(n -1) 4/5 (125)
However, this correlation should not be used in any calculation considering its confused origin.

Chen [61] performed an analytical study of laminar film condensation over a vertical bank of tubes.

He took into account the additional condensation between tubes and the effect of heat capacity:

h n h 1 = n -1/4 [
Though often mentioned, this work is never used for condenser design. Chen also noticed that the tem-530 perature difference could modify the inundation factor. He concluded that the smaller the temperature difference, the more important was the inundation. Asbik et al. [62] later adapted Chen's developments.

More recently, Murase et al. [63] obtained experimental results shown in Figure 23. They are close to Kern correlation and the present authors propose the following correlation from these results:

h n h 1 = n -1/7 (127)
and for the n-th tube:

535 h n h 1 = n 6/7 -(n -1) 6/7 (128)
In 2012, Ma et al. [64] developed an experimental procedure to artificially obtain the inundation factor of a vertical bank of tubes using an actually smaller one. Though no correlation was produced, it appears in Figure 24 that the experimental results lie above Kern line, which means it is slightly conservative.

Correlations based on mass rate 540

If previous correlations are based on the tube row number, which corresponds to the position of the tube within the bundle, the following ones are based on the condensate mass rate received by each tube.

The latter are therefore based on experimental results and are expressed as:

h n h 1 = Γ n γ n -s (129)
where Γ n is the condensation mass rate of the first n tubes, γ n is the condensation mass rate of the the n-th tube, and not the condensate mass rate draining onto the n-th tube. This way, the inundation factor equals 1 for n = 1. Γ n is therefore defined as:

Γ n = n i=1 γ i (130)
From experimental measurements performed on a rectangular bundle of 72 tubes divided into 11 rows with steam condensation, Fuks [65] obtained a value of s = 0.07. According to Bontemps [66],

Fuks did not manage to separate shear stress effects from inundation effects, therefore resulting in a low 550 inundation factor.

In 1968, Grant & Osment [67] performed an experimentation on an oval bundle of 139 staggered tubes with p/D o = 1.5 with low pressure steam. Using Fuks formulation, they obtained a coefficient s = 0.223. This formulation is commonly found in condensation literature.

The formulation of Wilson [68] is sometimes mentioned with a coefficient s = 0.16, which was 555 obtained by calculation to fit experimental data.

More recently, Hu & Zhang [69] have proposed a new inundation correlation with a variable coefficient s ranging from 0 at the top of the tube bundle to 0.37 at the bottom of it. According to the authors, the nil value at the top ensures an inundation factor equal to 1, which is frivolous since the ratio Γ 1 /γ 1 already equals 1. Regarding the value 0.37, it was obtain from the fitting of experimental results. the average inundation factor for each column, and then to calculate the mean value over the columns weighted by the number of tubes in each column.

Overall synthesis

In the case of industrial condensers, tube bundles may contains from hundreds to several thousands of 575 tubes which leads to an important condensate flow. Besides, considering that each tube acts as a sink for the vapour phase, the more tubes, the higher the vapour mass flow rate within the tube bundle.

Such a vapour flow is most likely to disturb the condensate as shown in Figure 26, making it modelling quite impossible.

If the vapour shear stress and the condensate inundation phenomena are both well enough described 580 when not occuring concomitantly, current literature is of little help when they do. The most frequent method consists of choosing a heat transfer coefficient for a single tube, and then to multiply it by an inundation factor.

Regarding the heat transfer coefficient on a single tube, the most appropriate correlation advised by the present authors is the one from Rose [52]. Indeed, it takes into account several phenomena without 585 being too conservative. However, the designer should be careful with this correlation, since it has been obtained for P ≤ 10. Were P be out of this range, the correlation of Fujii et al. [1] should be preferred.

As for the inundation factor, it appears that both the Kern equation [57] and Grant & Osment [67] are often used in condenser design. They give similar results, really close to experimental results being only slightly conservative. One is based on the tube row number while the other is based on the 590 condensate mass flow rate.

However, if the purpose of the modelling is to safely provide a thermal design of the condenser, one should use more conservative correlations for both the heat transfer coefficient and the inundation factor.

Once the local phenomena occurring near a tube in the bundle are acquired, the designer then needs 595 to figure out the tube bundle layout and how it will affect the vapour flow. This is a whole other matter that is rarely discussed by the industrials, since it is a valuable knowledge. for ( 12), ( 56) and ( 59)

2 . 3 .

 23 The vapour is saturated.85 The wall temperature T w is uniform.4. The condensate flow around the tube surface is laminar and purely orthoradial.5. The condensate film is not subject to shear stress at the liquid/vapour interface.6. The heat is transferred through the condensate film only by conduction. 7. The condensate properties are taken at saturation 90 In order to understand the further developments undertaken by reseachers after Nusselt's paper, it seems necessary to demonstrate his analytical work. Starting with the equation of momentum with the assumptions of a 2D stationary purely orthogonal flow within a slender film, one obtain the equation of movement for the condensate film:

  [41], who experimentally observed such a distribution. Then Memory & Rose solved Nusselt equation for various values of A in the range [0; 1]

270Figure 9

 9 Figure 9 describes the evolution of N u/ Re 1/2

285 20 Shekriladze&

 20 the sets of equation with both analytical and numerical tools. The analytical part shed light on new dimensionless numbers: the ρµ-ratio R and the condensation number H. Then, by means of Runge-Kutta-Gill method and for different values of R, H and F r, the equations were solved.Then they correlated the obtained results for high vapour flow. For F -→ 0, a dependence of N u to RH was highlighted as shown in Figure10. For large values of RH, the Nusselt number approach 290 Gomelauri solution (53) (-• line), while for small values of RH, the authors correlated the following Nusselt number (--line):

  ) are compared in Figure 11 for different values of RH. They both approach Nusselt solution for F -→ ∞ for all values of RH. Moreover, the formulation of Fujii et al. coincides with the formulation of Shekriladze & Gomelauri for RH -→ ∞.

  455

  ) and 20(b) for inline arrangements. When considering a staggered arrangement, tubes should be numbered as shown in Figure20(c).

5605. 3 Synthesis

 3 Using equations (105) and (106), the present authors have managed to compare both kinds of formulation. An iterative procedure was performed to calculate the mass rate of a vertical bank of tubes, in order to obtain the inundation factor of the n-th tube. Results are presented in Figure25with the inundation factors of Grant & Osment and Fuks. 565 The sharp decrease is also noticeable in the Grant & Osment dashed line. It is of great importance to notice that both the Kern formulation and the Grant & Osment formulation give similar results, since they are often recommended in condenser design literature. Therefore, we would recommend to use either of these two formulations : equations (111) and (112) for Kern or (129) with s = 0.223 for Grant & Osment.

  in steam condensers. Numerical Heat Transfer, Part A -Applications, 54(1):34-46, 2008.

Figure 1 :

 1 Figure 1: Different kinds of vapour and condensate flow in a shell-and-tube condenser

Figure 2 :Figure 3 :Figure 4 :

 234 Figure 2: Surface mass rate balance over a portion of condensate

Figure 5 :Figure 6 :Figure 7 :

 567 Figure 5: Crossflow within the condensate film

Figure 8 :Figure 9 :

 89 Figure 8: Separation point and recirculation flow

Figure 10 :Figure 11 : 47 Figure 12 :Figure 14 :Figure 15 : 50 Figure 16 :Figure 17 : 51 Figure 18 :Figure 19 :Figure 20 :Figure 21 :Figure 22 :Figure 23 :Figure 24 : 54 Figure 25 :Figure 26 :

 1011471214155016175118192021222324542526 Figure 10: Average Nusselt number for large oncoming vapour velocity

Table 1 :

 1 Comparison of the 3 expressions of the ratio N u/N u cp=0 for different Ja numbers These three expressions are compared in Table

			N u/N u cp=0
		1st approx. 2nd approx. Simplified
	Ja	(20)	(21)	(22)
	0.01	1.0020	1.0020	1.0020
	0.1	1.0197	1.0198	1.0198
	0.5	1.093	1.095	1.095
	1.0	1.175	1.180	1.183
	2.0	1.31	1.33	1.34
	3.0	1.43	1.45	1.48
	5.0	1.61	1.64	1.73

Table 2 :

 2 1/4 . Dependence of mean Nusselt number on a forced convection film condensation

	A N u Re	-1/2
	0.0	0.900
	0.2	0.906
	0.4	0.924
	0.6	0.953
	0.8	0.992
	1.0	1.040

Table 3 :

 3 Tube spacing parameter F d definition

	5	(122)

d + (1 -0.5647 F d ) n 4/5 -(n -1) 4/5 (123) According to equation (121), pure lateral drainage should be encountered for p t /D o ≤ 1.25 and pure 515 vertical drainage for p t /D o ≥ 1.40. In between, a combination of both inundation modes should be observed.

Uniform heat flux

Instead of assuming a cosine wall temperature, another approach was chosen by Fujii et al. [42]. The authors used the same assumptions as Nusselt, except the uniform wall temperature replaced by a uniform surface heat flux q. They obtained the following set of equations:

Using the same boundary conditions, the same mean velocity u is obtained. Therefore, substituting

For this experiment, boundary-layer separation occured at 83 • . So behind the separation point, the authors have neglected the friction force, solving the same equation as (1):

Then ( 44) and ( 45) are solved numerically by Runge-Kutta's method. Unfortunately, neither a local nor a mean Nusselt number was expressed from these results, only diagrams.

This method seems to neglect the fact that the vapour drag has an undesirable effect on condensation Therefore, this correlation should not be used for F > 10 in the current formulation. An asymptotic 405 model could be considered to correct it.

Vapour flow inclination effect

In 1974, Honda & Fujii [START_REF] Honda | Effect of the direction of on-coming vapor on laminar filmwise condensation on a horizontal cylinder[END_REF] studied analytically the case of inclined vapour flow within a plane orthogonal to the axis of the tube. As described in Figure 17, the oncoming vapour velocity direction has a angle φ with respect to the vertical. Based on the asymptotic shear stress of Shekriladze & Gomelauri [46],

the authors considered the following set of equations:

with the boundary-layer conditions:

= τ δ (shear stress at vapour-liquid interface)

Then they numerically solved the obtained differential equation for several values of F . Finally, they concluded that the average Nusselt number is slightly affected by the vapour flow orientation. Only the case of upward oncoming vapour velocity (φ > 5π/6) and F about unity has a significant effect on the 415 Nusselt number, since the condensate tend to flood the lower part of the tube.

Surface temperature variation

Carrying on with their previous work (see 3.4), Memory et al. [55] applied the same method in order to observe the effect of temperature distribution in the case of flowing vapour. For large oncoming velocities, their approach is based on Shekriladze & Gomelauri method [46]. They solved the differential 420 equation for several values of A (see (35)) and obtained the mean Nusselt numbers listed in Table 2.

For A = 0 (i.e. isothermal case), the original constant of Shekriladze & Gomelauri is obtained, but for increasing values of A, the heat transfer coefficient is slightly improved.

Influence of turbulence in the vapor and liquid phases. Journal of Heat Transfer, 121:874-885, 11 1999.