Solid-fluid-thermal analysis of heat exchangers


Brian Spalding

Keynote Lecture at 2005 ASME Heat Transfer Conference

July 17-22, 2005, San Francisco, California



Traditional heat-exchanger design methods predict steady-state uniform-property performance imperfectly; and
transient, varying-property flow, heat transfer, chemical reaction and stresses in solids not at all.

Traditional CFD techniques, relying on body-fitting grids and dubious turbulence models, handle only small-scale phenomena such as few-tube sub-sections of tube banks.

Intermediate CFD methods, which replace tube-bank details by space-averaged representations and heat-transfer correlations, have been successfully used for years; but they are little favoured by designers.


One reason is that heat-transfer correlations and other formulae needed by designers are too numerous and varied to be supplied with commercial codes; so users should provide Fortran themselves; but they lack time or skill.

What has long been needed is software which will itself understand formulae.

The lecture exemplifies how such software has been applied to a shell-and-tube heat exchanger.

It also shows how the stresses in the solid components of the heat exchanger can be computed at the same time.


Contents of the lecture

  1. Introduction
    1. Early successes
    2. Obstacles to further progress
    3. A way forward
  2. The Requirements of a Heat-Exchanger-Design Method
    1. Geometrical input data
    2. Material-property data
    3. Thermal and mass-flow boundary conditions
    4. Empirical correlations
    5. Predicting the flow pattern and temperature distribution
  3. Three ways of satisfying the requirements
    1. Method 1: 'User-supplied sub-routines'
    2. Method 2: 'Automated sub-routine writing'
    3. Method 3: When the program understands formulae

  4. Practical Example
    1. The heat exchanger and computer code in question
    2. Distributions of velocity, temperature, pressure and related properties
    3. Distributions of Reynolds, Prandtl and Nusselt number
    4. Variations of heat-transfer coefficient
    5. How the input data were introduced
    6. Discussion
  5. Simultaneous calculation of stresses in solids
    1. A common misconception
    2. The truth of the matter
    3. Fluid-structure interaction examples
    4. Exemplification for the heat exchanger
  6. Concluding remarks
  7. Acknowledgement
  8. References


1. Introduction

1.1 Early successes

The space-averaged-CFD technique applied to heat exchangers was first described than thirty years ago by Patankar and Spalding [1]. They concluded: "It therefore seems that a tool of considerable practical utility is in embryonic existence".

At first their expectations appeared to be fulfilled: application to two-phase flows helped to resolve difficulties encountered by the then growing nuclear-power industry.

Specifically, the shell-side steam-water mixture circulating in boilers, heated by pressurized water from the nuclear reactor, caused tubes to vibrate and baffles to corrode. Consequently, first US and German companies finally EPRI, sponsored the development of flow-simulating computer programs.

The work was pioneering and did not always proceed as rapidly as desired.

So one joker suggested that the name adopted for the EPRI-sponsored code, URSULA, was an acronym for:

 Urgently Required Solution Unacceptably Late Arriving.

Despite the implied criticism, the work was successful; and it was followed by the development of further computer codes for simulating steam condensers and cooling towers.

1.2 Obstacles to further progress

Nevertheless, the heat-exchanger-design community has shown little enthusiasm; and a recent paper [2] on the subject concluded "very few applications can be found of using CFD technique as a tool for heat-exchanger design optimization".

Instead, designers still prefer to use methods, for example those of Tinker [3] or Bell [4], in which the flow patterns are guessed , not calculated.

The reasons are not entirely clear. However, that they are in part psychological is suggested by the remarks of J Taborek [5] in the Hemisphere Handbook of Heat Exchanger Design.

He there opines: "Only if calculations are performed manually will the engineer develop a 'feel' for the design process as compared to the impersonal 'black box' calculations of a computer program".


1.3 A way forward

The approach recommended here may be congenial; for it enables designers to insert the same formulae, including Tinker-Bell 'correction factors', as they would supply to the 'hand-held calculators' preferred by Taborek.

The method to be described can be applied to heat exchangers of all types, to any participating fluids, and to any conditions of operation. However, in order to focus on essentials, discussion will henceforth be limited to:

Such a heat exchanger is shown from the outside here and in a sectioned view, but without the tubes, here. The software package allows the various components to be moved, re-sized, and viewed from various angles.

2. The Requirements of a Heat-Exchanger Design Method

2.1 Geometrical input data

Whatever the calculation technique, the apparatus in question must first be described in geometrical terms, including (for the simplest cases): next

2.2 Material-property data

Specification must be made of: However, for most materials, these properties are known to vary with temperature; and this knowledge is expressed either by way of:

If graphs are in question, their content must be converted into formulae or tables before it can be communicated to a computer program.

However, even when this has been done, the problem of using the information remains; for the whole point of a heat exchanger is to change temperature; and it is not known in advance what temperatures will prevail at any chosen point within the tubes or shell.

Therefore either:

The second is the main theme of the present paper.

2.3 Thermal and mass-flow boundary conditions

Also needed, of course, are the (known): The main task of performance prediction is to determine what will be the (mass-flow-weighted average) temperatures of the shell- and tube-side fluids at their outlets from the heat exchanger.

2.4 Empirical correlations

If: a traditional CFD code could predict performance well. However: This entails that, if performance were to be predicted purely from CFD, a very fine grid would have to be needed (see Prof. Sunden's lecture); and the performance-prediction calculation would take more time than any designer could afford to wait.

Moreover, for flow in tube banks, the knowledge of turbulence is still rudimentary; so the reliability of the predictions far from one hundred per cent.

The only practical solution is therefore to introduce additional information, derived from such experimental data as can be found, concerning the rates of heat and momentum transfer per unit area of solid-fluid interface.

This information is usually expressed in the form of mathematically-expressed relationships between the well-known 'dimensionless parameters'

namely: All these parameters involve the material properties listed in section 2.1.

Therefore the use of empirical correlations still requires the computer program to work out the property values from the given formulae and the temperatures which it finds at every point.

2.5 Predicting the flow pattern and temperature distribution

The outlet temperatures of the heat-exchanging fluids are the main quantities which it is desired to predict; however, these depend on the temperatures just upstream of the outlet; these just-upstream-of-outlet temperatures depend on the temperatures upstream of them; and so on.

Therefore, the whole temperature distribution has to be computed.

When the flow pattern is not one-dimensional, what point lies 'just upstream of' a given point is not obvious a priori;

therefore ability to calculate the temperature distribution depends on ability to calculate the whole three-dimensional flow distribution giving rise to it.


This therefore is what the computer program must additionally do, providing incidentally further information needed by the designer: the pressure losses suffered by the two streams.

Fortunately, computer programs (the so-called CFD codes) do exist for computing both the flow fields and the temperature distributions simultaneously.

Although their accuracy depends on the fineness of computational grid which is employed, and desirably fine grids do increase computer times and therefore costs, the requirements relating to shell-and-tube heat exchangers are usually affordable; but only when the space-averaged approach is adopted.


However, the heat-transfer and friction correlations require:

It follows that, even if the temperature variations were small enough not to affect material properties, the need for the code to evaluate formulae from values which varied from place to place would remain.


Thus the Reynolds Number enters most pressure-drop and convective-heat-transfer formulae; and its value depends on the local velocity, which varies with position in ways that are not known at the start.

Moreover, the angle between the mean-flow direction and the tube axis is probably important. It should appear in the formulae.

In summary, predicting the performance of shell-and-tube heat exchangers necessitates use of a program which is not restricted, as conventional design methods are, to the presumptions that:


3. Three Ways of Satisfying the Requirements

3.1 Method: 'User-supplied sub-routines'

Of course, many CFD codes have some built-in correlation-evaluation sequences, representing friction and heat-transfer; and they also contain computer-coding modules which express the variation with temperature and pressure of some properties of some materials.

In principle, there is no limit to the extent to which these provisions can be extended. But in practice, however much is provided, some users of the code will:


From the earliest years of commercial CFD, therefore, developers have allowed users to add coding modules of their own, usually in the form of Fortran or C subroutines, which would supplement the built-in correlations in the desired respect.

Users of the 1981 PHOENICS code, for example, will remember what clever use some users made of the so-called 'GROUND-coding' facility, which indeed many old-stagers continue to use. The paper by Stevanovich et al [2] is an excellent example of the use of this technique for heat exchangers.

However, the proportion of CFD-code users with the necessary skills is constantly diminishing; and few heat-exchanger designers either possess or have the time to acquire them.


3.2 Method 2: 'Automated sub-routine writing'

In order to enable PHOENICS users to benefit from the features of 'GROUND-coding' without themselves having to be familiar with either Fortran or C, the so-called 'PLANT' feature was introduced in 1997.

This enabled the user to express his wishes by way of formulae, written in accordance with prescribed rules; whereupon PHOENICS itself:


This was a big step forward; and it did, at least potentially, satisfy the 'formula-processing' requirement which has been pointed out above.

However, because it was not adequately presented to them, it did not convert many heat-exchanger designers into CFD users.

Perhaps also the 'prescribed rules' were shaped by those thinking too much of the Fortran to be written, and not sufficiently of the prospective user.

PLANT remains as a valuable feature of the current PHOENICS; but the feature now to be described has rendered it almost redundant.


3.3 Method 3: When the program understands formulae

The new feature has the acronym "In-Form", standing for "Input (or Intake) of Formulae". Its concept is very simple, namely: next

The formulae may provide CFD inputs of any kind, whether: The formulae most relevant to heat-exchanger design are those which: next

In-Form has many more capabilities than the heat-exchanger designer is likely to need. Therefore only the few relevant ones will be illustrated here, as follows: next

The relevance of these statements to heat-exchanger design should be obvious; for they say (nearly) all that is to be said about how local heat transfer from tube to fluid is to be calculated.

Once they have been typed into the input file, there is nothing further for the designer to do than to wait (for a minute or two) until the calculation has been completed.


4. Practical Example

4.1 The heat exchanger and computer code in question

A heat exchanger will be considered which has:

The shell-side fluid is water and the tube-side fluid SAE 5W-30 engine-oil ; these enter at 10 and 60 degrees Celsius respectively.

The mass flow rates of both fluids are 100 kg/s.

Text-book formula, connecting Nusselt, Reynolds and Prandtl numbers, have been adopted, with modifications made by the user, for the shell-side and tube-side heat transfer coefficients.

The computer code

The distributions of pressure, velocity and temperature of the two fluids have been computed by the general-purpose CFD code, PHOENICS, without activation of any heat-exchanger-specific special features.

Therefore any other CFD code equipped with a formula-understanding module could have been used.


The results to be displayed

The computer runs of which the results will be presented have been chosen in order to illustrate that:

The distributions are of course three-dimensional. However it is sufficient for present purposes to consider the central plane only.

Results will be represented by way of vector and contour plots on this plane.

After presentation of the results, some of the formulae used for generating them will be shown and discussed.


4.2 Distributions of velocity, temperature, pressure and related properties

The velocity distribution in the shell, on the central plane is shown in Fig. 1. The shell-fluid inlet is at the bottom on the left, and its outlet at the top on the right. The influences of the two baffles is clearly seen.

The temperature distributions in the shell- and tube-side fluids are shown in Fig.2 and Fig. 3. The latter conforms to that to be expected in an ideal counter-flow heat-exchanger, but the former (understandably) does not.

The pressure distribution in the shell is shown in Fig.4. It shows the to-be-expected discontinuities caused by the baffles.

As has been emphasized above, the properties of most fluids depend upon temperature. The variation of the viscosity of the tube-side fluid, for example, is shown in Fig. 5. Evidently (see scale on the right), it varies by a factor of more than two.

4.3 Distributions of Reynolds, Prandtl and Nusselt number

The dimensionless parameters which feature in the calculation of heat-transfer coefficients vary through the heat exchanger partly, because of temperature variations and partly because of the non-uniformities of velocity.

The distributions of Reynolds, Prandtl and Nusselt number calculated in the present case are shown, for the shell- and tube-side fluids in turn, in Fig. 6, Fig. 7, Fig. 8, Fig. 9, Fig. 10 and Fig. 11.

In respect of Reynolds number, variations are greatest on the shell side; the presence of the baffles accounts for this.

The consequential variation of Nusselt number is correspondingly great, being highest near the edges of the baffles and fifty times smaller in low-velocity regions, near the baffle roots and in the corners opposite inlet and outlet.

4.4 Variations of heat-transfer coefficient

From the Nusselt numbers, the computation of the corresponding shell- and tube-side heat-transfer coefficients per unit volume of shell is an easy step; and the overall heat-transfer coefficient, U, follows as their harmonic mean.

Their distributions are displayed in Fig.12, Fig.13 and Fig.14.

The shell-side coefficients show wide variations but, their effect on U is diminished by the fact that the tube-side ones are the smaller; so the maximum U is less than 50 % greater than the minimum U.

This observation might be thought to blunt somwhat the above criticism of conventional design methods; but it does not truly do so.

The objection to conventional methods is not that that they use a constant value of U for the whole heat exchanger; it is that the means of determining that value are based on guesswork and not on detailed analysis of the flow field.

4.5 How the input data were introduced

The results displayed so far will have come as no surprise to any heat-exchanger designer who has thought even superficially about what is truly happening inside his equipment.

What may interest him more may be to learn:

  1. how easily they were obtained, and
  2. how quickly the effects of changes to the input can be investigated.
The ability of the computer program to understand easily-written formulae is the key requirement.

Examples of such formulae will now be presented, starting with the following extract from the engineer's data-input file:


The lower-case words ('stored', 'var', 'is') tell the code that it must compute some new variables everywhere within the 'virtual heat exchanger' with which it is concerned. next

The upper-case words have been chosen by the engineer, who appears to have found:

The engineer has also chosen the formulae, most of which are rather easy to interpret. Thus:



How were the material properties calculated? Also by way of formulae drawn from textbooks. Thus, the formulae for the properties of the engine-oil in question were found in a book by K.Hagen [6] and translated thus into the form understood by (at least) the PHOENICS computer code: wherein POL4(TEM2,.....) signifies a fourth-power polynomial, the coefficients in which are the numbers which follow the comma.

Fortunately, many formulae are provided as items to be selected from a library. Nevertheless, it is not difficult to translate a new formula into 'In-Form-speak': users are not confined to what is on the code-vendor's menu.

4.6 Discussion

There are engineers who are unwilling to do more, when designing equipment, than make selections, by way of mouse-clicks, from what so many have done before that the code provider has had time and inducement to adopt as built-in options.

Others however, more like those referred to in the above quote from J.Taborek, will welcome the power which In-Form provides to "develop a 'feel' for the design process'.

Noting the considerable temperature differences between the metal temperature and the tube-fluid temperature, a designer might wonder whether he ought to multiply the relevant coefficient, as some textbooks recommend, by the ratio (wall viscosity/ bulk viscosity) raised to the power 0.14).

He could set his doubts at rest in a few moments by introducing the following lines into his input file
(by copying; God forbid that he should have to type it again!):

(stored var ENUM is: 10.^(POL4(TEMM,58.2987,-.53817,1.92827e-3, -3.16448,E-6,1.97922E-9)-2)

so as to compute the viscosity of the engine-oil at the wall temperature, TEMM; then he should modify the formula for coe2 so as to reflect the viscosity-ratio effect, thus:

(stored var COE2 is (ENUM/ENU2)^0.14*AOVERV*NUS2*CON2/DIAM).

If he did this, he would learn, after inspecting the computed results, that the total heat transferred would increase by only 0.5 %, so adding to his stock of 'engineering feel'.

Alternatively, he might wish to know how the heat-transfer rate per unit volume varied through the shell; then he could satisfy his curiosity by writing the single line:

(stored var flux is COEU*(TEM2-TEM))

and, one minute later, view the results presented in Fig. 15.

Comparing this with other images which he has seen before he will be able to decide whether he is satisfied with the nearly-two-fold variation there disclosed, or will seek by re-design to achieve greater uniformity.


Lastly, he might suddenly think: "What about fouling?" Then he has only to insert one line in his input file such as:

(stored var FOUL is .... whatever function of position and temperature he invents)

and to change the overall-coefficient line to:

(stored var COEU is 1/(1/COE1+1/COE2+FOUL))

The subsequent computer run will inform him of the consequences, both by way of numbers and of graphical displays.

There is no limit to the easily-effected variations of input and output which the In-Form facility allows; therefore, the main points having already been made, no further examples will be presented.


5. Simultaneous calculation of stresses in solids

5.1 A common misconception

It is widely believed that the problem of computing the stresses in solids differs from that of calculating fluid flow so much that different methods are needed for each problem.

In particular, it is often supposed that the so-called finite-element method must be used for the solid-stress problem.

Both these beliefs are erroneous.


5.2 The truth of the matter

In fact, the differential equations governing displacements in solids are very similar to (but simpler than) those governing velocities in fluids.

It is true that the displacements in the three coordinate directions are connected, via Poisson's Ratio, in a manner of which velocities know nothing.

Further, equations for the three components of rotation must be solved in addition

However, the finite-volume methods employed in CFD codes are well able to solve these equations, without any of the leftovers of pre-computer days (e.g. energy theorems and Galerkin weighting procedures) which clutter the finite-element formulation.

One of the first to recognise this was Stephen Beale [Reference 8].

5.3 Fluid-structure interaction examples

Results from an early study
are shown below for a two-solid-material block, heated by radiation from above, and cooled by a stream of air:

(1) velocity vectors,

(2) displacement vectors, computed at the same time, and

(3) horizontal-direction stresses, obtained by post-processing.


A bent beam

A beam bent by fluid forces is shown here.

Flow in a turn-around duct
One of the early arguments for preferring FE methods was that they handled curved boundaries more easily.

But FV techniqes serve just as well, even when the grid is cartesian, as shown here for flow in a turn-around duct.

In this case there are heat sources in the solids, giving rise to these thermal-expansion contours.

Because of the unified-computational-mechanics technique, we are able to compute simultaneously the velocity and displacement vectors.

Forces on an under-water structure

Finally, deflection is shown of a sea-bed structure, resulting from the action of
surface waves.

Here it is essential to be able to simulate solids and fluids at the same time.


5.4 Exemplification for the heat exchanger

In order to illustrate these assertions, calculations have been performed for the longitudinal stresses in the tubes of the present shell-and-tube heat exchanger.

Fig. 16 shows the results of the calculations in terms of vector diagrams.

Vectors of what?

  1. The left half of the diagram shows vectors of shell-side-fluid velocity; and
  2. the right-half shows vectors of displacement in the tube metal.
The latter have been shifted to the right of the former for ease of display. However PHOENICS specialists might recognise that the MUSES (Multiply-USed-Space) trick has been played, whereby a double-size grid has been employed, with:

Why do the displacement vectors have the distributions which have been shown? Because of the arbitrarily-introduced fixings postulated for the tube ends within the shell. The stressescan be deduced from the strains, which are themselves deduced by differentiating the displacements.

Fig. 17 shows contours of the stresses, of course in the right-hand part of the display region.

The higher (compressive) stresses appear in the lower part of the shell, where the tube ends have been rigidly fixed.

The importance of Fig. 17 in the present context is that:


6. Concluding remarks

It should be stressed that the techniques employed can be applied also: next

It is therefore argued that:
  1. the time has come for heat-exchanger designers to exploit computational fluid dynamics with space-averaged representations of heat-transfer and friction processes;

  2. they should consider calculating solid stresses at the same time;

  3. specifically, the difficulty of introducing the necessary empirical knowledge has been removed by the formula-parsing capability;

  4. this enables the computer code to 'understand' the formulae which designers are accustomed to use in their traditional methods;

  5. thus all that formerly could be done only by specialists capable of exploiting the 'user-subroutine' facility , can now be done by any flexibly-minded heat-transfer engineer.

7. Acknowledgement

The author gratefully acknowledges the assistance of Mr. Nikolay Pavitsky, of CHAM's office in Moscow, Russia, in developing the In-Form facility during the last five years.


8. References

  1. SV Patankar and DB Spalding "A calculation procedure for the transient and steady behavior of shell and tube heat exchangers" in "Heat exchangers: design and theory sourcebook" edited by N Afgan and E Schluender, McGraw Hill 1974
  2. Z.Stevanovic, G Ilic, N Radojkovic, M Vukic, V Stefanovic and G Vukovic "Design of shell and tube heat exchangers by using CFD technique - part one: thermo-hydraulic calculation". Facta Universitatis Series Mechanical Engineering vol 1 no 8, 2001, pp1091 -1105
  3. T Tinker J. Heat Transfer vol 80 pp 36-52 1958
  4. KJ Bell "Final report of the cooperative research program on shell-and-tube heat exchangers" University of Delaware Exp.Sta.Bull. 5 1993
  5. J Taborek "Recommended method: principles and limitations" in "Hemisphere Handbook of Heat Exchanger Design" ed. by GF Hewitt, Hemisphere, New York 1983
  6. K Hagen "Heat transfer, with applications", Prentice Hall, 1999
  7. D Chisholm. "Two-phase flow in shell-and-tube heat exchangers" in "Heat Exchanger Technology" ed. by D.Chisholm, Elsevier, 1988
  8. S B Beale and SR Elias "Numerical solution of Two-Dimensional Elasticity problems by means of a 'SIMPLE-based' finite-difference scheme". NRC Report No 32090, 1991/02

    The End !!!