Transfer function-based analysis of the frequency-domain properties of a double pipe heat exchanger

This paper discusses irrational transfer function representation for a typical thick-walled double pipe heat exchanger. Closed-form analytical expressions for the individual elements of 2 × 2 transfer function matrix are derived both for the parallel- and the counter-flow configurations using the Laplace transform method. Based on the transfer function representation, its frequency responses are demonstrated both in the form of three-dimensional graphs as well as classical, two-dimensional Bode and Nyquist plots. Finally, steady-state temperature profiles are presented and compared for both flow arrangements considered.


Introduction
The transfer of thermal energy between fluids is one of the most important processes in engineering, which is possible through the use of different types of heat exchangers. They are broadly used in many thermal processes, both for the heating and cooling operations, and one of the simplest are the so-called double pipe heat exchangers. Despite their simplicity and low efficiency, they are very important due to their educational value in teaching the basics of heat exchanger design. Moreover, they form the structural basis for many more complex constructions, such as e.g. shell and tube heat exchangers or steam generators [24].
The steady-state properties of this class of thermal devices are very important from the point of view of their operation and control and thus are readily available since decades in many professional textbooks and scientific publications. However, since heat exchangers usually constitute part of larger industrial systems, their transient operations can occur frequently and thus can significantly influence overall system performance [21,22,25]. Therefore, contemporary mathematical models of heat exchangers should comprise not only their steady-state behaviour but also their dynamic properties, which are essential e.g. to impart effective control and to take preventive measures considering the safety aspect of the entire plant [5,11].
The conventional tool for the mathematical description of the so-called distributed parameter systems (DPSs), which also include heat exchangers, are partial differential equations (PDEs) [13]. The development of the PDE-based mathematical modeling methods has progressed here in two parallel directions, with the use of either numerical or analytical approaches. For example, the numerical approach to some problems of plate fin and tube heat exchangers is presented in [23], whereas the same collection also includes the closed-form analytical solutions for the transient heat and moisture diffusion in a double-layer plate [8]. Moreover, there exist many professional software packages which allow to obtain numerical solutions to PDEs and make it possible to solve even very complex mathematical models, including the non-linear ones [19].
However, the use of such applications is usually limited to the personal computers or workstations. They are usually of little use for the real-time applications such as hardwarein-the-loop, or for implementation on embedded systems. From the control theory point of view, this kind of analytical representation is not fully effective because it does not express directly the relationships between the output (controlled) and the input (manipulated) variables of the system, and thus it is difficult to be used e.g. in the control system synthesis. In contrast, analytical methods allow for more thorough physical interpretation of the process and are more useful in the automatic control theory [20].
In this paper a method of the analytical description of a double pipe heat exchanger is considered, which is based on its transfer function representation. The main advantage of this approach, as compared to the numerical modelling, is that it constitutes a convenient starting point for computer implementation of different control algorithms. The knowledge of the transfer function make it to possible to determine not only the steady-state but also frequency-and time-domain responses of the system. As opposed to the rational transfer function describing the dynamic properties of lumped parameter systems (LPSs), transfer functions of DPSs have the form of irrational functions [7,27]. Due to the mathematical complexity, their analysis is more difficult, and possible applications are more limited than in the case of the finitedimensional models. Therefore, in order to enable the implementation of the developed over the years techniques for the synthesis of control systems, the infinite-dimensional DPS models are usually replaced by their finitedimensional approximations [3,4,9,12,16]. Nevertheless, regardless of the approximation method used, the starting point for the synthesis of a control system should be based on the possibly most accurate description of the DPS, taking into account its infinite-dimensional nature, e.g. a model in the form of the irrational transfer function.
Compared to the previously presented results, in the present paper the thermal capacity of the internal tube of the exchanger is taken into account, and two different typical flow configurations are considered, both for the steady-state and frequency domains. The remainder of the paper is organized as follows. Section 2 introduces a mathematical model of the considered thick-walled double pipe heat exchanger in the form of hyperbolic PDE system. In Sect. 3, transfer function matrices for the parallel-and counter-flow configurations are derived and analyzed. Based on the obtained transfer functions, Sect. 4 presents selected frequency responses of the exchanger, both in the form of three-dimensional graphs as well as classical, two-dimensional Bode and Nyquist plots. Next, Sect. 5 compares the steady-state temperature distributions for both flow configurations considered. Finally, short conclusions and future work prospects are given in Sect. 6.

Governing PDEs
The considered heat exchanger consists of two concentric pipes (tubes) containing fluids flowing from the inlet of each tube towards its outlet (Fig. 1). In order to avoid ambiguity, throughout the rest of the paper the external tube will be referred as shell and the internal-simply as tube. Heat is transferred from one fluid to the other through the wall of the tube, either from tube side to shell side or vice versa. Depending on the flow arrangement, the fluids can enter the shell and the tube from the same or from the opposite ends of the exchanger. The first configuration is commonly known as parallel-while the second is usually referred to as counter-flow one.
In order to develop the mathematical model of the exchanger based on the energy balance equations, the following simplifying assumptions are made: • exchanger is perfectly insulated from the environment; • there are no internal thermal energy sources; • the flows are sufficiently turbulent to cause effective heat transfer; • only forced heat convection is considered (i.e. longitudinal heat conduction within the fluids and wall is neglected); • pressure drops of fluids along the shell and the tube are negligible; • the densities and heat capacities of the shell, tube and fluids are time and space invariant; • the convective heat transfer coefficients are constant and uniform over each surface.
According to the above assumptions, double pipe heat exchanger depicted in Fig. 1 is governed, based on the thermal energy balance equations, by the following PDE system [2,6,18]: where t 2 ½ 0; þ1Þ represents time, l 2 ½ 0; L stands for the space variable and the constant parameters k 1 , k 2 , k 3 , k 4 depend on the shell and tube diameters, physical parameters of the fluids and the exchanger material [18]: where q is the density, c is the specific heat and h is the heat transfer coefficient (subscripts as for the temperatures, see Fig. 1). As one can notice, Eqs. (1-3) are weakly coupled, i.e. coupled only through the terms of convective heat exchange which do not contain derivatives. In contrast to the commonly adopted heat exchanger models such as those considered e.g. in [1,11,17], the abovementioned simplifications do not contain any assumption on the negligible thermal capacity of the tube, which results in the additional Eq. (3) representing heat conduction through the wall of the tube.
In order to obtain a unique solution of Eqs. (1-3), one must specify the appropriate initial and boundary conditions. The initial conditions can be specified in the following form: where # t0 ðlÞ; # s0 ðlÞ; # w0 ðlÞ : ½0; L ! R are given functions representing the initial (i.e. determined for t ¼ 0) temperature profiles along the heat exchanger.
The form of the boundary conditions for Eqs. (1-3) depends on the flow arrangement (see Fig. 1). For the case of the parallel-flow one obtains the following boundary conditions: whereas for the counter-flow configuration the boundary conditions are expressed as follows: where # ti ðtÞ, # si ðlÞ represent the inlet temperatures which can be considered, from the control theory point of view, either as control signals or external disturbances.

Transfer function representation
For the heat exchanger under consideration one can distinguish two lumped input signals represented by the boundary conditions (6) or (7), and two output signals representing the distribution of the fluid temperatures along its axis. As shown e.g. in [7,10], it is possible for the considered DPS with boundary inputs to obtain closed form expressions for the transfer functions by taking the Laplace transforms of the original PDEs and solving the resulting boundary value problem. In the following two subsections, transfer function matrices for the heat exchanger are introduced for two different boundary input configurations corresponding to the parallel-and counter-flow modes. Next, closed-form analytical expressions for those transfer functions are derived.

Parallel-flow mode
The transfer function matrix of the heat exchanger can be defined for the parallel-flow configuration described by the boundary conditions (6) for # t ð0; sÞ ¼ 0, all for zero initial conditions, stands for the Laplace transform of #ðl; tÞ in variable t.
The vector of Laplace transforms of the fluid temperatures can therefore be determined, assuming zero initial conditions, based on the following equation: where can be considered as the vector of the boundary input signals for the parallel-flow mode.
A block diagram of the transfer function model for this mode is presented in Fig. 2. In order to find the elements of the transfer function matrix Gðl; sÞ in (8), the Laplace transform approach can be used to solve Eqs. (1-3) with the boundary conditions (6). After Laplace transformation (11), Eqs. (1-3) take the following form: Assuming zero initial conditions and expressing the Laplace-transformed wall temperature # w ðl; sÞ from (17) as the equation set (15)(16)(17) can be reduced to the following two equations:  where p tt ðsÞ ¼ À p ss ðsÞ ¼ À can be considered as elements of the following matrix: PðsÞ ¼ p tt ðsÞ p ts ðsÞ p st ðsÞ p ss ðsÞ By performing the Laplace transform again, now with respect to the space variable l: and taking into account that one can transform the Eqs. (19) and (20) Finding the inverse Laplace transform of (28) and (29) with respect to q by taking advantage of the following property [14]: where SðqÞ and TðqÞ represent polynomials in q of degree M and N [ M, respectively, and k j is a single root of TðqÞ, yields the following form of the equations: with p tt ðsÞ, p ts ðsÞ, p st ðsÞ and p ss ðsÞ being the elements of the matrix PðsÞ in (25) given by Eqs. (21)(22)(23)(24) and / 1 ðsÞ, / 2 ðsÞ-its eigenvalues given by (31-33).

The case of zero heat transfer
As previously mentioned, the hyperbolic PDEs (1-3) are weakly coupled through the terms on the right-hand sides which represent the convective heat exchange. Assuming that no heat is transferred between the fluids through the , one obtains the extremely simplified form of the equations for which the matrix PðsÞ in (25) takes the following form: with the eigenvalues which consequently result in the following transfer function matrix: for the parallel-flow, and for the counter-flow configuration. The resulting system can thus be considered as two separate pure time-delay subsystems with time delays given by (59) or (60), representing two fluids of constant temperature profiles travelling along the exchanger.

Frequency responses
Based on the transfer functions derived in Sect. 3 it is straightforward to determine frequency responses of the heat exchanger. For this purpose, one should replace in the relationships (37-40) or (49-52) the operator variable s with the expression ix, where i is the imaginary unit and x stands for the angular frequency.
The graphical representation of these responses can take the form of three-dimensional graphs, taking into account their dependence on both the angular frequency x and the space variable l. Another possibility is the representation in the form of classical two-dimensional plots, determined for fixed value of the space variable, e.g. l ¼ 0 or l ¼ L.
Considering as an example the Bode plot, the expressions for the logarithmic gain and phase take the following wellknown form [14]: and u l; where the expressions for the modulus and argument of the frequency response are as follows: From the practical point of view, most important are the responses evaluated at the exchanger outlets, i.e. assuming l ¼ 0 for " G st ; " G ss , and l ¼ L for the remaining transfer functions. Fig. 6 shows classical two-dimensional Bode plots of the heat exchanger frequency responses G ts ðl; sÞ determined for l ¼ L. Next, the Nyquist plot for the same transfer function channel is presented in Fig. 7.
As seen from the Bode plots, the increase in the frequency of the sinusoidal input signal initially causes a decrease in the amplitude of the output signal, and then it gives rise to a local maximum. To a lesser extent this phenomenon is apparent in the phase characteristic. As the frequency increases, this effect repeats itself, which can be observed as oscillations on the Bode plot or as characteristic ''loops'' on the Nyquist plot. These oscillations are closely associated with the wave phenomena taking place inside the exchanger pipes [14].
Analysis of the frequency responses of the heat exchanger exhibits typical characteristics of systems with distributed delay [15]. In particular, in the case of the  ''straightforward'' transfer function channels G tt ðl; sÞ, G ss ðl; sÞ and their counter-flow counterparts, one can notice the dominant influence of the transport delay in the fluid flow. In order to confirm this conclusion, Nyquist plot of the frequency response G ss ðL; sÞ is shown in Fig. 8. The amplitude damping of the sinusoidal oscillations in the real and imaginary parts of the frequency response is relatively small, which is reflected here in the circular-shaped plot.
On the other hand, in the case of the ''crossover'' transfer functions G ts ðl; sÞ and G st ðl; sÞ, the damping of the input signal with increasing frequency is much greater as for the ''straight-forward'' channels (see Figs. 6, 7).

Steady-state analysis
The steady-state temperature profiles along the heat exchanger can be determined directly from Eqs. (1-3) by assuming all time derivatives equal to zero, and solving the resulting boundary value problem with the boundary conditions representing constant inlet temperatures for the parallel-flow or for the counter-flow configuration The other possibility is to calculate the steady state responses of the exchanger from its transfer functions assuming s ¼ 0 or, equivalently, from its frequency responses assuming x ¼ 0, and this approach will be applied below. For s ¼ 0 one obtains the following values of the parameters (21)(22)(23)(24): p ss ¼ p ss ð0Þ ¼ À for which one obtains based on (32) and (33) and consequently from (31) Based on the transfer functions determined in Sect. 3, it is now straightforward to obtain the formulas for the steadystate temperature profiles:   Figure 9 shows the steady-state temperature profiles for the tube-and shell-side fluids, calculated based on Eqs. (78-83) for the following constant values of the inlet temperatures: From the obtained results it is possible to determine e.g. the outlet temperatures of both fluids involved in the heat exchange. For example, for the parallel-flow configuration the outlet temperature # t ðLÞ of the heated fluid is about 54:8 C and the outlet temperature # s ðLÞ of the heating fluid is 69 C (see Fig. 9). Reducing the flow rate v t of the heated fluid from 1 to 0:2 m/s increases its outlet temperature # t ðLÞ to about 68:6 C and also causes an increase in the outlet temperature # s ðLÞ of the heating fluid to 75:9 C. As is apparent from Fig. 10, the change in the flow configuration causes further increase in the temperature of the heated fluid as compared to the parallel-flow mode. For example, when changing the flow rate v t from 1 to 0:2 m/s its outlet temperature # t ðLÞ reaches 71:3 C.
It is worth noting that the derived analytical expressions describing the steady-state properties not only make it possible to determine the outlet temperatures of the fluids but also allow an analysis of the temperature profiles along the heat exchanger, which may be of great importance from a technological point of view.
To sum up, the counter-flow mode has several advantages as compared to the parallel-flow one. The outlet temperature of the heated fluid can approach the inlet temperature of the heating fluid. The more uniform temperature difference between the two fluids prevents thermal stresses in the exchanger material. The other advantage is that the more uniform difference between temperatures # t ðlÞ and # s ðlÞ has an effect of more uniform heat transfer rate. The results presented above have been compared and found to be in general consistent with those well-known from the literature [18,24,26].

Conclusion
The closed-form analytical expressions have been derived for the individual elements of the 2 9 2 transfer function matrix of a double pipe heat exchanger working both in the parallel-and counter-flow modes. Unlike the case of lumped parameter systems, the transfer functions derived for this distributed parameter system contain irrational functions such as exponential and square root ones. As shown in the paper, the location of the boundary inputs of the system related to the flow configuration significantly affects its transfer function representation. Based on the obtained transfer functions, selected frequency responses of the heat exchanger have been presented both in the form of three-dimensional graphs as well as classical, twodimensional Bode and Nyquist plots. Moreover, steady state temperature profiles have been determined and compared for both considered flow arrangements. The future works could include determination of the space-time impulse responses for the individual transfer function channels, both for the parallel-and the counterflow mode. Another issue to be thoroughly examined, which is very important from the control synthesis point of view, is selecting an appropriate method for the transfer functions and space-time responses approximation using the finite-dimensional models mentioned in Sect. 1.
Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and the source are credited. Fig. 9 The steady state temperature profiles # t ðlÞ, # s ðlÞ and # w ðlÞ for the parallel-flow heat exchanger (# si ¼ 100 C, # ti ¼ 50 C) Fig. 10 The steady state temperature profiles # t ðlÞ, # s ðlÞ and # w ðlÞ for the countercurrent-flow heat exchanger (# si ¼ 100 C, # ti ¼ 50 C)