Access the full text.

Sign up today, get DeepDyve free for 14 days.

Mathematical and Computer Modelling of Dynamical Systems
, Volume 14 (3): 18 – Jun 1, 2008

/lp/taylor-francis/mathematical-modelling-of-a-domestic-heating-system-with-stratified-HBMcxIXtl6

References for this paper are not available at this time. We will be adding them shortly, thank you for your patience.

- Publisher
- Taylor & Francis
- Copyright
- Copyright Taylor & Francis Group, LLC
- ISSN
- 1744-5051
- eISSN
- 1387-3954
- DOI
- 10.1080/13873950701844907
- Publisher site
- See Article on Publisher Site

Mathematical and Computer Modelling of Dynamical Systems Vol. 14, No. 3, June 2008, 231–248 Mathematical modelling of a domestic heating system with stratiﬁed storage tank a a b T. Kreuzinger *, M. Bitzer and W. Marquardt Mechatronic Engineering, Corporate Sector Research and Advance Engineering of Robert Bosch GmbH, Robert – Bosch – Strasse 2, 71701, Schwieberdingen, Germany; Lehrstuhl fu¨r Prozesstechnik, RWTH Aachen University, Turmstraße 46, 52056, Aachen, Germany (Received 30 October 2006; ﬁnal version received 17 August 2007) A hybrid distributed parameter model of a heating system for domestic hot water is presented in this paper. This heating system comprises a condensing boiler (burner), a counter current heat exchanger, and a so-called stratiﬁed storage tank which is the state of the art domestic hot water storage unit. The paper presents the model for the diﬀerent operational modes of the plant which are described by a ﬁnite state automaton representing the discrete-event dynamics and driving the underlying continuous-time dynamics of the storage tank, the heat exchanger, and the burner. These interconnected components are modelled by a system of six coupled, quasi-linear partial diﬀerential equations (PDEs) comprising diﬀusion-, convection-, and source terms. In order to perform numerical simulations, the set of PDEs is spatially discretized using the method of lines. Thereby, the inﬂuence of various discretization schemes on the temporal evolution of the traveling temperature proﬁles in the single components is investigated. A high resolution slope limiter scheme for the stratiﬁed storage tank and a higher order up–/downwind scheme for the heat exchanger and the burner are found to be an appropriate choice for the spatial discretization of the model equations in order to adequately cover the plant dynamics. Simulation results fortify the eﬀectiveness of the chosen discretization schemes and show the excellent performance of the suggested model representing the measurement data. Keywords: stratiﬁed storage tank; hybrid distributed parameter system; diﬀusion- convection system; counter current heat exchanger 1. Introduction The use of storage tanks for domestic hot water in contemporary heating systems, as shown in Figure 1, has proven great potential to save energy and increase comfort and availability. For solar heating systems, it is even indispensable, since in most cases there is a signiﬁcant time gap between energy gain and consumption. Although storage of hot water is common practice in most households, the eﬃciency of ordinary storage systems is not satisfactory. Therefore, in the heating systems community, the state-of-the-art system for domestic hot water storage is the so-called stratiﬁed storage tank [1]. In general, the ﬂuid in the storage vessel has a certain temperature proﬁle and thus a density gradient *Corresponding author. Email: tobias.kreuzinger@de.bosch.com th Revised and expanded version of a paper presented at the 5 Vienna International Conference on th Mathematical Modelling (5 Mathmod), Vienna, Austria, 2006. ISSN 1387-3954 print/ISSN 1744-5051 online 2008 Taylor & Francis DOI: 10.1080/13873950701844907 http://www.informaworld.com 232 T. Kreuzinger et al. which assures the physical principle of stratiﬁcation. To maintain this stratiﬁcation, the storage tank in Figure 2 is driven with relatively low mass ﬂows and equipped with baﬄe plates to ensure that the mixing and excitation of agitation ﬂows in the in- and outﬂow regions of the tank are minimized. Due to the characteristic stratiﬁcation of the water in the tank, i.e. hot water at the top part and cold water at the bottom part, domestic hot water (e.g. for shower, bath, kitchen) is always taken from the uppermost ﬂuid layer. During tapping, equal amounts of cold ﬂuid are introduced at the very bottom of the tank such that the vessel is always entirely ﬁlled with water. During loading, the storage tank is connected to the secondary side of a counter current heat exchanger and ﬁlled with the desired amount of hot water. Thereby, cold water is pumped from the bottom part over the heat exchanger and introduced back up into the top part of the tank. To feed the corresponding primary side with hot water, the heat exchanger is connected to a condensing boiler, referred to as the burner. The interconnection of the components of Figure 1. Heating unit with stratiﬁed storage tank, control unit, loading pump, burner, counter current heat exchanger, and temperature sensors (NTC 1 and NTC 2). Figure 2. Sketch of the plant. Mathematical and Computer Modelling of Dynamical Systems 233 interest, i.e. burner, heat exchanger, and storage tank, can be inferred from Figure 2. The system considered here is referred to as the Cerasmart Module and commercially available from Robert Bosch GmbH. The dynamics of the plant is characterized by discrete-time events due to e.g. tapping or loading events as well as continuous time-varying spatial temperature proﬁles, which occur in the storage tank, the heat exchanger, and the burner. These components are represented by the respective distributed parameter models which are then coupled according to the ﬂowsheet in Figure 2. In detail, the models are given as a diﬀusion– convection system (DCS) with heat loss to the ambiance for the storage tank and as two counter current transport systems with heat exchange between the primary and secondary side for the heat exchanger and the burner. Thus, the entire model consists of a ﬁnite state automaton interacting with a continuous-time model of six underlying coupled quasi- linear partial diﬀerential equations. Consequently, the entire plant can be interpreted as a hybrid distributed parameter system [2]. Due to the rapid increase of the cost of energy, consumers require that process control schemes operate the plant in such a manner that energy consumption is minimized. However, today’s standard is that the loading cycle is driven by a two-point controller that activates the loading pump and the burner if a predeﬁned temperature threshold at the temperature sensor NTC1 is detected and correspondingly stops the loading after detecting a respective temperature threshold at NTC2. To develop a process control system which minimizes the consumed energy of the heating system, the set-up of a model-based process control scheme [3] is mandatory. Even though complex 3D ﬁnite element models exist for similar applications [4], the contribution of this paper is to provide a reduced-order model of the plant which is suitable for the development of a process control scheme for energy eﬃcient, optimal loading of the storage tank. The tasks therefore are the design of a state observer [5], a control algorithm based e.g. on a ﬂatness approach [6,7], or on the method of characteristics [8,9] and a dynamic real-time optimization scheme for the design of an economically optimal loading strategy [10]. The ﬁrst part of the paper outlines the derivation of the hybrid plant model. Here, the ﬁnite state automaton with its ﬁve discrete states and its corresponding transitions is presented. Furthermore, the partial diﬀerential equations (PDEs) for the storage tank, the heat exchanger and the burner are derived. To perform numerical simulations, the spatial discretization of the PDEs is investigated in the second part of the paper. Thereby, the temporal evolution of the traveling temperature proﬁles in dependence of the discretization approach for the convection terms is analysed. For the tank model, a ﬁnite volume approach switching between an up– and a downwind scheme in dependence of the ﬂow direction is compared to a switching high resolution slope limiter scheme. For the heat exchanger and the burner, a ﬁnite volume ﬁrst order up- and downwind approach is compared to respective third- and fourth-order discretization approaches. In the last part of the paper, the suggested plant model is identiﬁed and compared to measurement data. 2. Hybrid plant model The discrete plant dynamics is characterized by switching events due to tapping and loading, which drive the underlying dynamical behaviour of the heat exchanger, the burner, and the storage tank [11]. The discrete-event dynamics is fully described by a ﬁnite state automaton of ﬁve discrete states and their respective transitions. These states represent the 234 T. Kreuzinger et al. operational modes of the plant and are distinguished as follows: idle state, loading state, tapping state, and the combined state of tapping and loading. Thereby, the operational mode where tapping and loading occurs simultaneously has to be divided into two discrete states depending on the resulting direction of the plug ﬂow in the tank. The underlying continuous-time dynamics is characterized by time-varying spatial temperature proﬁles in the storage tank, the heat exchanger, and the burner. These components are modelled by six quasi-linear PDEs for the six state variables, namely the temperature T(z,t) in the storage tank, the primary and secondary side temperatures Y (x,t)and Y (x,t) in the heat pr sec exchanger, and the primary side, wall, and secondary side temperatures t (z, t), t (z, t), pr w and t (z, t) in the burner. These temperature proﬁles depend on the spatial coordinates z, sec x, and z as well as on time t. In the subsequent section (2.1), the ﬁnite state automaton is introduced, followed by a description of the distributed energy balances for the heat exchanger, the burner, and the storage tank in Sections 2.2–2.4. 2.1. Discrete state automaton A ﬁnite state automaton can be denoted as a labelled transition system, described by a triple (L,A,E), where L is the ﬁnite state space, A is a ﬁnite set called alphabet, and E is the corresponding transition rule [2]. A sequence (l , a , l , a ,..., l , a , l )with(l , a , 0 0 1 1 N71 N71 N m m l ) 2 E for m ¼ 1, 2, . . . is called a trajectory on the state space [2]. Using this notation, mþ1 the state space of the ﬁnite state automaton of the storage tank is given by L ¼fl ¼ idle; l ¼ tap; l ¼ load; l ¼ load > tap; l ¼ load < tapg: ð1Þ 1 2 3 4 5 The graph of the resulting ﬁnite state automaton is depicted in Figure 3. Note that depending on the tapping and loading mass ﬂows during the states l and l , the resulting 4 5 plug ﬂow in the storage tank is either positive, i.e. directed from top to bottom, or negative, i.e. from bottom to top. For an appropriate choice of the numerical approach, the separation of these two states is of great relevance. The alphabet of the ﬁnite state Figure 3. Graph of the ﬁnite state automaton describing the operational modes of the plant. For l ;l m n the sake of simplicity, transitions are labelled only unidirectionally, i.e. for each t there exists a l ;l m n corresponding transition t . k Mathematical and Computer Modelling of Dynamical Systems 235 l ;l m n automaton [2] is labelled by implicitly caused discrete switching times t (due to a tapping or loading event) denoting the transitions from one state l to another state l , i.e. m n no l ;l m n A ¼ t ; with l ; l 2 L; k ¼ 1; 2; ... ð2Þ m n l ;l m n With Equations (1) and (2), the general transition rule can be derived to be ðl ; t ; l Þ m n l ;l l ;l l ;l 1 2 2 4 4 3 [2], i.e. ðl ; t ; l ; t ; l ; t ; l ; ...Þ is an example for a possible trajectory on the discrete 1 2 4 3 1 2 3 state space L in dependence of the events occurring. 2.2. Distributed parameter model of the heat exchanger The counter current heat exchanger is modelled by a set of two coupled ﬁrst-order partial diﬀerential equations [8,12], each comprising a convection and a source term accounting for the heat exchange between the primary and the secondary side of the system. Therefore, a pipe-in-pipe system with length L is used to model the primary (outer pipe) and the secondary side (inner pipe) of the heat exchanger (see Figure 4). The separating wall between the two pipes is assumed to be negligibly thin, such that heat conduction along and through the wall is neglected. Setting up the energy balance for an inﬁnitesimally small volume element with length dx and applying a Taylor series expansion for the element bound at x þ dx yields the model equations of the heat exchanger for dx ! 0. The resulting equation for the primary side is @Y ðx; tÞ @Y ðx; tÞ 2ar pr pr sec ¼ v ðtÞ Y ðx; tÞ Y ðx; tÞ ; h pr sec 2 2 @t @x rðYÞcðYÞ r r pr sec x 2½ð 0; L Þ; t < t < t: 3Þ Accordingly, the energy balance of the secondary side is given by @Y ðx; tÞ @Y ðx; tÞ 2a sec sec ¼v ðtÞ þ Y ðx; tÞ Y ðx; tÞ ; l pr sec @t @x rðYÞcðYÞr sec x 2ðð 0; L ; t < t < t: 4Þ Figure 4. Sketch of the counter current heat exchanger with primary and secondary side temperature proﬁles Y (x,t) and Y (x,t) as well as the respective ﬂow velocities v (t) and v (t). pr sec h l 236 T. Kreuzinger et al. As shown in Section 4, the heat transfer coeﬃcient a between the two heat exchanger pipes in Equations (3) and (4) is obtained by minimizing the least squares error between model and measurement data. Therefore, the geometry of the pipes is chosen to match the transfer area and the primary and secondary side volumes of the heat exchanger. The Dirichlet boundary condition (BC) at x ¼ 0 reads as in Y ð0; tÞ¼ Y ðtÞ; t < t < t ð5Þ sec sec and analogously at x ¼ L as in Y ðL ; tÞ¼ Y ðtÞ; t < t < t: ð6Þ pr Y pr The initial conditions (ICs) in the pipes are 0 0 Y ðx; tÞ¼ Y ðxÞ;Y ðx; tÞ¼ Y ðxÞ; x 2½0; L : ð7Þ pr sec Y pr sec l ;l l ;l m n m n Here, t ¼ t ; m 2f1; 2g; n 2f3; 4; 5g; and t ¼ t ; m 2f3; 4; 5g; n 2f1; 2g, which k k means that the heat exchanger model is activated during the discrete states l , l 3 4 and l . 2.3. Distributed parameter model of the burner This section ﬁrst summarizes the operating principle of the condensing boiler, followed by the presentation of the distributed parameter model of the burner. A fan feeds the combustion chamber of the burner with an air mass ﬂow proportional to the fan speed (n (t)) (see Figure 5). The respective gas mass ﬂow is adjusted by a mechanical l – control, which guarantees optimal combustion in the burner. Here, the combustion chamber is referred to as the primary side of the burner. The burning gas is then used to heat up the water in the secondary side of the burner. The model of the burner is therefore set up in accordance to the counter current principle used for the heat exchanger in Section 2.2. Due to the thickness of the walls in the burner, i.e. the separating wall between the primary and secondary side, the structure in Equations (3) and (4) has to be extended by an additional PDE for the wall. As depicted in the sketch in Figure 5, the actual model of the burner is split up into two parts. The ﬁrst part consists of two algebraic relations to determine the Figure 5. Setup of the burner composed of the fan with the fan speed n and the heat exchanging unit (pipe system) with the primary and secondary side temperatures t (z, t) and t (z, t) as well as pr sec the respective ﬂow velocities v (t) and v (t). g h Mathematical and Computer Modelling of Dynamical Systems 237 in ﬂow velocity v (t) and the primary side feed-ﬂow temperature t ðtÞ as a function of the fan pr speed n (t) according to n ðtÞ v ðtÞ¼ v ; ð8Þ g g;max and in t ðtÞ¼ þ p : ð9Þ p 3 pr 4 n ðtÞ þ p F 2 The parameters v and p are determined by minimizing the least squares error g,max 1,...,4 between simulations and measurement data, i.e. all eﬀects of combustion are summarized in in t ðtÞ. Equations (8) and (9) are employed to calculate the respective primary side inputs pr of the second part of the burner model, i.e. the dynamical part (the pipe system) which can be interpreted as a counter current heat exchanger with a dividing wall, as depicted in Figure 6. Setting up the energy balance for the volume element dz and following the procedure in Section 2.2 yields the model equations for the primary side, @t ðz; tÞ @t ðz; tÞ pr pr ¼ v ðtÞ g t ðz; tÞ t ðz; tÞ ; z 2 ½0; L Þ; t < t < t; ð10Þ g pr w t pr @t @z the wall, @t ðz; tÞ g w pr sec ¼ t ðz; tÞ t ðz; tÞ ðÞ t ðz; tÞ t ðz; tÞ ; pr w w sec @t q q pr sec z 2½0; L ; t < t < t; ð11Þ and the secondary side, @t ðz; tÞ @t ðz; tÞ sec sec ¼v ðtÞ þ gðÞ t ðz; tÞ t ðz; tÞ ; z 2 ð0; L ; t < t < t: ð12Þ h w sec t sec @t @z Figure 6. Sketch of the heat exchanging part of the burner with the primary side, wall, and secondary side temperatures t , t , and t as well as the respective ﬂow velocities v (t) and v (t). pr w sec g h 238 T. Kreuzinger et al. The boundary conditions at the inlets of the primary and secondary side are in in t ðL ; tÞ¼ t ðtÞ; t ð0; tÞ¼ t ðtÞ; t < t < t: ð13Þ pr t sec pr sec The ICs read as 0 0 t ðz; tÞ¼ t ðzÞ; t ðz; tÞ¼ t ðzÞ; z 2½0; L : ð14Þ pr sec t pr sec In the sequel, the IC of the wall t (z, t) is chosen to be the arithmetic mean between the temperatures of the primary and secondary side. Note that due to the signiﬁcant length of the pipe system in the burner, heat conduction within the wall is negligible. The coeﬃcients of the exchange terms in Equations (10)–(12), i.e. g , g , q and q , are also obtained by pr sec pr sec an optimization-based parameter identiﬁcation. The optimization result is presented in Section (4). 2.4. Distributed parameter model of the stratiﬁed storage tank The governing model equations of the storage tank can be obtained following the derivation of the heat exchanger and burner equations. Due to the large diameter of the storage vessel, ﬂow velocities are signiﬁcantly smaller than in the heat exchanger. Therefore, it is necessary to extend the partial diﬀerential equation by a conductive component which is set up according to Fourier’s heat conduction [13]. Further, the BCs have to be expanded to Cauchy BCs to model the heat loss through the ceiling and bottom of the storage tank. The resulting temperature distribution in the stratiﬁed storage tank can be described by the following PDE with BCs and a respective IC according to @Tðz; tÞ @ Tðz; tÞ @Tðz; tÞ PDE : rðTÞcðTÞ ¼ l ðTÞ rðTÞcðTÞvðtÞ @t @z @z 4k ðT; vÞ DTðz; tÞ; z 2ð0; L Þ; t < t < t ; ð15Þ @Tð0; tÞ in BC : 0 ¼ l ðTÞ rðTÞcðTÞvðtÞ Tð0; tÞ T ðtÞ @z v ðtÞ in þ rðTÞcðTÞ T ðtÞ Tð0; tÞ k DTð0; tÞ; t < t < t ; ð16Þ @TðL ; tÞ v ðtÞ T d in 0 ¼ l ðTÞ rðTÞcðTÞ T ðtÞ TðL ; tÞ w T @z a þ k DTðL ; tÞ; t < t < t ; ð17Þ b T IC : Tðz; t Þ¼ T ðzÞ; z 2½0; L ; ð18Þ with DT ¼ T – T , the density and heat capacity of water r(T)and c(T) as well as the amp loading and tapping velocities v (t)and v (t). Note that due to the high Peclet numbers of 1 d about 5 10 , the temperature dependency of r, c and lw is negligibly small in the Mathematical and Computer Modelling of Dynamical Systems 239 considered range of operation. The respective starting and end times are deﬁned as l ;l l ;l m n m n t ¼ t ; m 2f1g; n 2f2; 3; 4; 5g; and t ¼ t ; m 2f2; 3; 4; 5g; n 2f1g. Equations (15)– k k (18) hold for the discrete states l – l . Further note that due to the occurrence of complex 2 5 free convection phenomena in the state l , this case is not analysed further in this contribution. The notation for the loading and tapping velocities and the resulting plug ﬂow velocity v(t) ¼ (v(t)– v (t))/a can be inferred from Figure 7, where the parameter a describes the ratio between the cross-sectional areas of the storage tank and the heat exchanger pipes. For the presented model, the parameter l (T) only covers the heat conduction due to intermolecular transport, i.e. l ¼ l . Note, however, that w w,mol due to the forced convection, macroscopic back-mixing phenomena are present at the storage wall and contribute to the heat conduction coeﬃcient such that the actual l 4 l [14]. w w,mol The general hybrid I/O-behaviour of the plant is on the one hand dominated by uncertain tapping events driven by the consumer and on the other hand based on the loading events initiated by a respective process control strategy. These events determine the in- and outﬂows of the storage tank. Consequently, the model equations depend on the corresponding discrete states in L and may furthermore be simpliﬁed signiﬁcantly. This concerns especially the BCs (16)–(17) of the stratiﬁed storage tank. The reduction for the relevant discrete states in L is summarized below and can be easily incorporated in the respective equations according to l : v ðtÞ¼ 0; v ðtÞ¼ av ðtÞ; and vðtÞ < 0; 2 1 d l : v ðtÞ¼ av ðtÞ; v ðtÞ¼ 0; and vðtÞ > 0; 3 1 d v ðtÞv ðtÞ l : v ðtÞ 6¼ 0; v ðtÞ 6¼ 0; and vðtÞ¼ > 0; 4 1 d v ðtÞv ðtÞ 1 d l : v ðtÞ 6¼ 0; v ðtÞ 6¼ 0; and vðtÞ¼ < 0: 5 1 d Figure 7. Left: labelled scheme of the storage tank with respective in- and outﬂow velocities v (t) l,d and the resulting plug ﬂow velocity v(t). Right: cross section of the storage wall with corresponding diameters d , m ¼ 1,2,3, heat conduction coeﬃcients l , and heat transfer coeﬃcients a , n2{i, a}. m m n 240 T. Kreuzinger et al. 2.4.1. Model parameters In the sequel, the derivation of the heat transfer coeﬃcient of the storage tank is outlined. The derivation of the respective coeﬃcient is hereby based on physical relations. As shown in Figure 7, the storage wall consists of three material layers with diﬀerent heat conductivity. The coeﬃcient k is therefore composed of the convective heat transfer in water and air, the heat conduction in each of the three wall layers and the radiative heat at the outside wall of the tank. According to [13], the equation for the heat transfer coeﬃcient k of a wall with various layers and an inner wall diameter d is given by "# 1 d d d mþ1 k ¼ þ ln þ ; ð19Þ a 2l d a d i m m a mþ1 m¼1 with d and l denoting the m-th diameter and the m-th conductivity, see Figure 7, right. m m Inside the tank, the heat transfer coeﬃcient reads as a (T,v) ¼ Nu l /d, with the i w 2 0,333 corresponding Nusselt number Nu(T,v) ¼ (49,028 þ 4,173Pr v(t) d /(nL )) for water, assuming laminar ﬂow and suﬃcient friction between the ﬂuid and inside wall of the tank [15,16]. The heat transfer coeﬃcient at the outside wall of the tank is composed of a convective and a radiant part yielding a ¼ a þ a . The convective part can a convection radiation be calculated by a ¼ Nu l /LT with the Nusselt number convection air 0;296 2 0;563 0;167 Nu ¼ 0; 825 þ 0; 387 Ra 1 þð0; 493=PrÞ ð20Þ at the ambient and the Rayleigh number [13] L gðr r Þ T amb air;wall Ra ¼ : ð21Þ na r air air;wall The radiative heat transfer is governed by 8 3 a ¼ 5; 67 10 E 4T ð22Þ radiation amb with the emission ratio E [15]. The derivation of these equations assumes that the outside wall temperature is approximately equal to the ambient temperature. This can easily be aﬃrmed for suﬃciently good isolation. For a profound discussion on the appropriate choice of parameters for diﬀerent materials and geometries, the reader is referred to [16]. 3. Numerical simulation of the hybrid plant model This section illustrates the numerical setup of the hybrid plant model for the discrete states l – l . Due to the hyperbolic nature of the PDEs during tapping and loading, the 2 5 discretization of the convection terms plays a major role for numerical simulation studies. Mathematical and Computer Modelling of Dynamical Systems 241 Hence, diﬀerent discretization schemes are investigated in terms of accuracy and resulting model dimension of the discretized lumped model. The ﬁrst part of this section considers the discretization of the storage tank model. Here, the reconstruction of sharp spatial transitions, as they are present in the stratiﬁed storage tank, is of great interest. Therefore, a ﬁnite volume up- and downwind scheme in dependence of the ﬂow direction as well as a high resolution scheme, a so-called slope limiter scheme [17], are investigated. The second part of the section gives a brief discussion on the discretization of the heat exchanger and burner model where up- and downwind schemes of diﬀerent approximation orders are used. 3.1. Discretization of the storage tank equation A method of lines (MOL) approach [18] is employed for the spatial discretization of the model Equation (15). The temperature distribution in z – direction is approximated by clustering the spatial domain in N ﬁnite volumes, see Figure 8. The centre of the i-th volume is denoted by z and the volume bounds as z 1 and z 1, respectively. Integration i iþ 2 2 of (15) over a control volume yields z z iþ1=2 iþ1=2 Z Z @T l @ T @T 4k w w dz ¼ vðtÞ DT dz: ð23Þ @t rc @z @z drc z z i1=2 i1=2 For the evaluation of the integral (23), certain assumptions are made. Choosing the coeﬃcients l , r, c, k and the temperature T to be constant within the interval z 2 [z w w i71/ , z ] and with dz ¼ z – z ¼ const, entails the discretized version of the storage 2 iþ1/2 iþ1/2 i71/2 PDE (15), i.e. i iþ1 i i1 iþ1=2 i1=2 dT l T 2T þ T ½T T 4k w w ¼ vðtÞ ðT T Þ: ð24Þ amb dt rc dz dz drc As mentioned before, special eﬀort has to be spent to ﬁnd an appropriate approximation of the convective terms. The various approaches mentioned previously are generally represented by the following substitution for the control volume bounds i+1/2 T according to i1=2 i1=2 i1 i i1 v > 0 : T ¼ T þ ðT T Þ; ð25Þ Figure 8. Finite volumes of discretized boundaries () and PDE (.). 242 T. Kreuzinger et al. iþ1=2 iþ1=2 i iþ1 i T ¼ T þ ðT T Þ; ð26Þ i1=2 i1=2 i i i1 v < 0 : T ¼ T ðT T Þ; ð27Þ iþ1=2 iþ1=2 iþ1 iþ1 i T ¼ T ðT T Þ; ð28Þ where F denotes the so-called limiter function [17,19]. Consequently, for v 4 0, the spatially discretized PDE (15) results as hi iþ1=2 i1=2 i F iþ1 i i1 F i i1 i iþ1 i i1 T þ ðT T Þ T ðT T Þ dT l T 2T þ T 2 2 ¼ vðtÞ dt rc dz dz 4k ðT T Þ; i ¼ 1::N : ð29Þ amb T drc In Equation (29), the well-known upwind scheme is obtained by setting the limiter iþ1/2 i71/2 function F ¼ F ¼ 0, whereas the corresponding downwind scheme can be iþ1/2 i71/2 obtained for F ¼ F ¼ 2 [18]. Since the plug ﬂow velocity in the storage tank changes sign, Gibb’s phenomenon [20] has to be avoided by employing the upwind scheme for v 4 0 and the downwind scheme for v 5 0, i.e. the choice of the numerical scheme depends on the ﬂow direction and therefore also on the discrete state of the ﬁnite state automaton. The so-called central diﬀerence scheme (CDS) is obtained by setting iþ1/2 i71/2 F ¼ F ¼ 1. The performance of the previously named approaches changes signiﬁcantly for diﬀerent spatial resolutions, especially for the upwind scheme. However, due to numerical diﬀusion, this approach is not feasible for the reconstruction of sharp spatial transitions. As opposed to the upwind scheme, the CDS perfectly reconstructs sharp transitions of the spatial proﬁle. Unfortunately, this discretization method leads to the occurrence of spatial oscillations and therefore to an unsatisfactory result. The idea of a high resolution slope limiter scheme is now to combine the beneﬁts of the previously named methods, i.e. oscillation-free reconstruction of sharp transitions. This can be realized by deﬁning F as a nonlinear function of the so-called smoothness y, which in turn is deﬁned as i1 i2 i1=2 T T y ¼ i i1 T T v > 0 : ð30Þ i i1 iþ1=2 T T y ¼ iþ1 i T T iþ1 i i1=2 ~ T T y ¼ i i1 T T v < 0 : ð31Þ iþ2 iþ1 iþ1=2 T T y ¼ iþ1 i T T for the corresponding ﬂow directions. In the literature, many limiter functions can be found [19,21]. Numerical tests were performed for a number of limiter functions, such as the monotonized centred, Van Leer, Koren, and Superbee limiters [17,19,21]. For the Mathematical and Computer Modelling of Dynamical Systems 243 purpose of reconstructing sharp spatial transitions, the Superbee limiter with the limiter function ~ ~ ~ FðyÞ¼ max 0; min ð1; 2yÞ; min ð2; yÞ ð32Þ yields the most promising results (see Section 4). It is furthermore worth noting that the presented discretization schemes are conservative [19]. 3.2. Discretization of the heat exchanger and burner equations Similar to the spatial discretization of the stratiﬁed storage tank, a MOL approach is used for the discretization of the heat exchanger and the burner. Due to the similarity of the system structure and the occurring spatial temperature proﬁles in the heat exchanger and the burner, the subsequent considerations are outlined on the basis of the heat exchanger model in Equations (3) and (4) but are then used for the discretization of both the heat exchanger and the burner. The following three discretization approaches adopted from [18] are investigated: j j1 @Y Y Y 1st order : ¼ ð33Þ @x dx jþ1 j j1 j2 @Y 2Y þ 3Y 6Y þ Y 3rd order : ¼ ð34Þ @x 6dx jþ1 j j1 j2 j3 @Y 3Y þ 10Y 18Y þ 6Y Y 4th order : ¼ : ð35Þ @x 12dx Hereby, a classical ﬁrst-order up/downwind scheme is used for the calculation of the boundary values in Equations (33) and (34), while the forth-order approach is based on a special higher order up- and downwind approximation of the boundary conditions (referred to as ‘‘Subroutine DSS020’’ in [18]). To classify the transient and stationary performance of the diﬀerent discretization orders, a respective reference trajectory has to be speciﬁed. This can be done using an analytical solution of the PDE. However, the analytical solution is generally not known or diﬃcult to evaluate. For the special case of the heat exchanger, the analytical solution includes inﬁnite sums of Bessel functions [22]. Their evaluation is costly and only possible to a certain accuracy. Hence, an alternative approach is suggested for cases, where an exact inversion of the PDE considered is available [6]. The procedure is then as follows: a desired trajectory w(t) which is dynamically feasible on one hand and which covers the transient and the stationary behaviour of the process on the other hand has to be chosen. With the use of an exact inversion S of the system S, a feedforward control u(t) can be calculated, i.e. u(t) ¼ S w(t) [15]. Then w(t) can be used as the respective reference trajectory when the calculated feedforward control is applied to the system S, i.e. y(t) ¼ Su(t) ¼ SS w(t) and therefore yðtÞ ¼ wðtÞ (see Figure 9). Substituting the distributed system S by its discretized counterpart yields the distorted output j¼1::N Y 244 T. Kreuzinger et al. P P trajectory y ~ðtÞ¼ wðtÞ, which is then used for the assessment of the j¼1::N discretization scheme. For the explicit case of the heat exchanger model, the performance is depicted in Figure 10. The inversion of the heat exchanger model is hereby obtained according to the methodology described in [6]. Figure 10 further shows that the minimum variation between the desired and distorted trajectory can be obtained using the 3rd order discretization scheme. This means that for the considered range of operation, the diﬀerent approximation of the boundary values in Equation (35) leads to an overall degradation of the 4th order scheme compared to the 1st and 3rd order approach. Hence, the 3rd order scheme according to Equation (34) is used for the spatial discretization of both, the heat exchanger and the burner model. 3.3. Simulation results of the hybrid plant model The discretized plant model, i.e. burner, heat exchanger, and stratiﬁed storage tank, is implemented in a MATLAB/SIMULINK environment [23], allowing for the numerical simulation of the plant and the examination of its hybrid behaviour. The implementation is set up in a ﬂexible manner, such that all types of scenarios which can possibly occur in the course of the plant operation may be studied. Figure 11 shows the simulation result of the temperature proﬁle T(z,t) in the storage tank for various variations of the ﬂow Figure 9. Structure for the classiﬁcation of discretization schemes using exact system inversion of the distributed parameter system S. Figure 10. Left: transient and stationary distortion due to spatial discretization of the heat exchanger PDE using N ¼ 50 compartments. Right: magniﬁed version showing the transition between the transient and stationary region. Mathematical and Computer Modelling of Dynamical Systems 245 velocities v and the fan speed n . Hereby, the simulation scenario covers the following l,d F eﬀects: at t ¼ 0 s the storage tank has an initial temperature proﬁle with T(z,0) ¼ 293 K. During the time interval 0 s t 5 700 s, the storage tank is operated in the discrete state l with a ramp-like rise of the fan speed n for 100 s t 300 s and a loading velocity 3 F v ¼ 0.3 m/s, according to the solid lines in Figure 11. It can be seen that due to the wall and the circular setup of the heat exchanger and the burner, the ramp-like increase of the in fan speed results in a smooth rise of the feed ﬂow temperature T ðtÞ at z ¼ 0. Due to the tapping event occurring during the time interval 700 s t 1100 s, cold water is introduced at z ¼ L . Note that in this case, the discrete state automaton is in state l , T 5 where the tapping ﬂow velocity v exceeds the loading ﬂow velocity v , i.e. the plug ﬂow d l reverses and the tank discharges. Hereby, the transitions between the discrete states l and l are denoted by the arrows in Figure 11. At t ¼ 1100 s, the state automaton switches back to the state l and the storage tank gets reﬁlled. For this latter case, the variation of the loading velocity for t 1070 s is chosen such that the variation of the feed ﬂow in temperature T ðtÞ at the top of the storage tank is kept inside a certain band. The respective region is marked by a circle in Figure 11. 4. Model identiﬁcation In this section, the hybrid plant model is identiﬁed via least squares optimization using measurement data obtained for the operational modes l and l . The datasets of the heat 2 3 exchanger, the burner, and the storage tank are presented and compared to the respective simulation results. Thereby, the discretization scheme and the required number of grid points for the storage tank model are furthermore discussed and determined. It shall be emphasized that the subsequent plots rely on a normalized time scale. Figure 12 shows the identiﬁcation results for the heat exchanger and the burner. It can be inferred from the plots that the heat exchanger and the burner model are feasible to track the measurement data very well for both, the stationary and the transient regions. To identify the storage tank model, ten equidistantly distributed thermocouples are attached from top to bottom on the inside wall of the tank. The numerical simulations of the upwind scheme Figure 11. Simulation results of the plant model for loading and tapping scenarios. Top left: tapping and loading velocities v and v . Bottom left: variation of the fan speed n . Right: 3D d l F temperature proﬁle T(z,t) of the storage tank. 246 T. Kreuzinger et al. (Figure 13) and the slope limiter scheme using a Superbee limiter (Figure 14) are compared to corresponding measurement data. In the case of the upwind scheme with a spatial discretization of 20 nodes, the model reconstructs the measured values very poorly. However, by increasing the spatial discretization up to 300 nodes, a far better ﬁt can be achieved. Even though better results are attained by increasing the spatial resolution of the upwind scheme, numerical diﬀusion is still present, see Figure 13. As illustrated in Figure 14, the slope limiter scheme bears signiﬁcant improvements for the loading scenario. It is furthermore shown that applying the high resolution scheme with Superbee limiter to the tapping state where hot water is taken from an entirely ﬁlled storage cell, the measured and simulated values agree very well. Note that for the high resolution scheme, a spatial discretization of only 21 grid nodes is suﬃcient. Figure 12. Left: plot of the measured and simulated outlet temperatures of the primary and secondary side of the counter current heat exchanger. Right: measured and simulated temperature proﬁles at the secondary outlet of the condensing boiler. Figure 13. Measured and simulated temperature proﬁles for a loading scenario, using an upwind diﬀerence scheme for discretization of the PDE with N ¼ 20 grid points (left) and N ¼ 300 grid T T points (right). Mathematical and Computer Modelling of Dynamical Systems 247 Figure 14. Plot of the measurements and simulations using the Superbee limiter and N ¼ 21 grid points for the discretization of the PDE. Left: loading scenario l . Right: tapping scenario l . 3 2 5. Conclusion and outlook In this paper, the derivation of a hybrid model of a domestic heating system comprising a condensing boiler (burner), a heat exchanger, and a stratiﬁed storage tank is outlined. Thereby, the model parameters of the burner and the heat exchanger are obtained by minimizing the least squares error between the simulation results and measurement data. The derivation of the storage tank parameters are based on physical considerations. To perform numerical simulations, diﬀerent discretization schemes according to the method of lines approach are investigated in terms of feasibility and numerical performance. It is hereby shown that for the storage tank, where sharp spatial temperature gradients are present, the high resolution slope limiter scheme with Superbee limiter outperformed the up-/downwind scheme and the CDS. It is further found that a third-order up-/downwind scheme is well suited for the discretization of the heat exchanger and the burner. The good agreement between the measurement data and simulation results for stationary and transient scenarios justiﬁes the suitability of the suggested models of the burner, the heat exchanger, and the stratiﬁed storage tank. Consequently, this model can be used to study the challenging task of advanced process control design for domestic heating systems with stratiﬁed storage tank as an example for a hybrid distributed parameter system [24]. Future research is necessary to model the free convection phenomena occurring in the course of the discrete state l of cooling. Note 1. The MATLAB - ode23tb solver is used for the time-integration of the discretized PDEs. References [1] I. Dincer, and M.A. Rosen, Thermal Energy Storage - Systems and Applications, John Wiley, New York, 2002. [2] A. van der Schaft, and H. Schumacher, An Introduction to Hybrid Dynamical Systems, Springer Verlag, Berlin, 2000. [3] W.H. Ray, Advanced Process Control, McGraw-Hill Education, New York, 1980. [4] W. Lin, and S.W. Armﬁeld, Direct simulation of natural convection cooling in a vertical circular cylinder, Int. J. Heat Mass Trans. 42 (1999), pp. 4117–4130. 248 T. Kreuzinger et al. [5] T. Kreuzinger, M. Bitzer, and W. Marquardt, Design of a quasi-linear distributed state observer and parameter estimator for a stratiﬁed storage tank, in Proceedings of the 2006 IEEE International Conference on Control Applications, Munich, Germany, 2006, pp. 680–685. [6] T. Kreuzinger, F.A. King, M. Bitzer, and W. Marquardt, Feedforward boundary control of linear distributed parameter systems using a numerical inverse Laplace transform, preprints of the 3rd IFAC Symposium on System, Structure and Control SSSC07, Foz do Iguassu, Brazil, 17–19 October 2007. [7] J. Rudolph, Randsteuerung von Wa¨rmetauschern mit o¨rtlich verteilten Parametern: ein ﬂachheitsbasierter Zugang, at-Automatisierungstechnik 48 (2000), pp. 399–406. [8] H.-K. Rhee, R. Aris, and N.R. Amundson, First-Order Partial Diﬀerential Equations, Vol. I, Prentice-Hall, Englewood Cliﬀs, NJ, USA, 1986. [9] Huilan Shang. Characteristic-based Control of Distributed Parameter Systems, Ph.D. thesis, Department of Chemical and Materials Engineering, University of Alberta, 2002. [10] J. Busch, J. Oldenburg, M. Santos, A. Cruse, and W. Marquardt, Dynamic predictive scheduling of operational strategies for continuous processes using mixed – logic dynamic optimization, Comput. Chem. Eng. 31 (2007), pp. 574–587. [11] T. Kreuzinger, M. Bitzer, and W. Marquardt, Hybrid modeling of a stratiﬁed storage tank,I. Troch, and F. Breitenecker, eds., in Proceedings of the 5th Vienna Symposium on Mathematical Modelling (5th Mathmod), Vienna, Austria, 2006, Vol. 2, pp. 1–10. [12] L.C. Evans, Partial Diﬀerential Equations, American Mathematical Society, Rhode Island, USA, 2002. [13] R.B. Bird, W.E. Stewart, and E.N. Lightfoot, Transport Phenomena, John Wiley, New York, USA, 1960. [14] R.F. Probstein, Physiochemical Hydrodynamics, Wiley Interscience, New York, 2003. [15] B. Gluck, Wa¨rmeubertragung - Wa¨rmeabgabe von Raumheizﬂa¨chen und Rohren. VEB Verlag fur ¨ ¨ ¨ Bauwesen, Berlin, Germany, 1989. [16] VDI Gesellschaft Verfahrenstechnik und Chemieingenieurwesen (ed.), VDI – Wa¨rmeatlas. VDI Verlag GmbH, Du¨ sseldorf, 1991. [17] P.L. Roe, Some contributions to the modelling of discontinuous ﬂows, Lecture Notes Appl. Math. 22 (1985), pp. 163–193. [18] W.E. Schiesser, The Numerical Method of Lines Integration of Partial Diﬀerential Equations, Academic Press, San Diego, 1991. [19] R. LeVeque, Numerical Methods for Conservation Laws, Birkha¨ user Verlag, Basel, 1992. [20] S.V. Patankar, Numerical Heat Transfer and Fluid Flow, McGraw-Hill, New York, 1980. [21] R. Ko¨ hler, Preprocessing Tool for Method-of-Lines Discretization of Process Models in the Simulation Environment DIVA, Fortschritt-Berichte Reihe 20 Nr. 357, VDI-Verlag Du¨ sseldorf, [22] M.A. Jawson and W. Smith, Countercurrent transfer processes in the non-steady state, in Proceedings of the Royal Society, vol. 226, 1954. [23] The MathWorks Inc. Matlab Documentation, August 2005. Matlab 7.1, Release 14. [24] P.D. Christoﬁdes, Control of nonlinear distributed process systems: recent developments and challenges, AIChE J. 47 (2001), pp. 514–518.

Mathematical and Computer Modelling of Dynamical Systems – Taylor & Francis

**Published: ** Jun 1, 2008

**Keywords: **stratified storage tank; hybrid distributed parameter system; diffusion-convection system; counter current heat exchanger

Loading...

You can share this free article with as many people as you like with the url below! We hope you enjoy this feature!

Read and print from thousands of top scholarly journals.

System error. Please try again!

Already have an account? Log in

Bookmark this article. You can see your Bookmarks on your DeepDyve Library.

To save an article, **log in** first, or **sign up** for a DeepDyve account if you don’t already have one.

Copy and paste the desired citation format or use the link below to download a file formatted for EndNote

Access the full text.

Sign up today, get DeepDyve free for 14 days.

All DeepDyve websites use cookies to improve your online experience. They were placed on your computer when you launched this website. You can change your cookie settings through your browser.