Real-time state observer design for solar thermal heating systems

Real-time state observer design for solar thermal heating systems

Applied Mathematics and Computation 218 (2012) 11558–11569 Contents lists available at SciVerse ScienceDirect Applied Mathematics and Computation jo...

971KB Sizes 0 Downloads 43 Views

Applied Mathematics and Computation 218 (2012) 11558–11569

Contents lists available at SciVerse ScienceDirect

Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc

Real-time state observer design for solar thermal heating systems R. Kicsiny ⇑, Z. Varga }, Hungary Institute of Mathematics and Informatics, Szent István University, Gödöllo

a r t i c l e

i n f o

Keywords: Real-time state observer Solar thermal heating system Linear mathematical model Estimation algorithm Simulation results Measurements

a b s t r a c t Nowadays it is important to investigate and develop solar thermal heating systems as an environmental friendly technology. In the paper we demonstrate a physically-based linear mathematical model that can be applied to a wide range of solar heating systems to analyze and develop them. We construct a Luenberger type state observer for the mathematical model, which enables us to estimate in real-time the unmeasured state variables of the solar heating system. We give a general algorithm for the practical application of the estimation process. Comparing calculations and measurements, we introduce the usability of the state observer in case of a particular solar heating system that produces domestic hot water at the campus of the Szent István University, Hungary. Ó 2012 Elsevier Inc. All rights reserved.

1. Introduction It often occurs in technical practice that all variables determining the state of some system are not measured (as a function of time), because for example the needed measurement cannot be carried out physically (see e.g. [1]), may be either too complicated, or too expensive (see e.g. [2]). In such a case the estimation of the unmeasured state variables is needed. In case of systems with either discrete or continuous time mathematical model, the problem can be solved in many cases by the methods of mathematical systems theory (for basics of this discipline see [3], for a more recent reference see [4]), by determining a so-called state observer [5]. If the estimation of the unmeasured state variables is performed in real-time, then we can not only monitor the system but also apply a feedback control strategy to it, e.g. for stabilization or optimization. In case of the so-called completely observable systems, the state of the system can be asymptotically reconstructed by the determination of an appropriate state observer. For continuous-time linear systems (that are described by systems of ordinary linear differential equations) different approaches are available, also depending on whether the problem is stochastic or deterministic. For an application of a state observer for a linear system in engineering (and for further references) see e.g. [1]. There the investigated stochastic mathematical model describes the dynamics of continuous fluidized bed dryer processes. The method used for estimating the state variables, in order to establish an optimal feedback control, is based on the Kalman filter, and the study proposes a procedure for the on-line estimation of the moisture content of the dried material. In [6] the discrete time mathematical model of a simple hydraulic system with two state variables and with random noise is introduced to illustrate the discussed mathematical method. One of the state variables is measured by which the other one is estimated with a state observer. For error filtering different types of the Kalman filter are used here, as well. Unlike the above studies, our model will be deterministic. ⇑ Corresponding author. E-mail addresses: [email protected] (R. Kicsiny), [email protected] (Z. Varga). 0096-3003/$ - see front matter Ó 2012 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.amc.2012.05.040

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

11559

State observers corresponding to linear time-varying systems are discussed in [7] and [8]. The first one discusses minimal order observers, the latter deals with full order observers that are applicable under particular conditions. Our dynamic model, instead, in accordance with the discrete character of the measurements, will be a piecewise autonomous linear system. In [9] a state observer is obtained for the time-invariant linearization of a nonlinear mathematical model of a ZnS precipitation process, based on the observability theory of linear systems. The following mathematical model, that will be later applied in our work as well, is related to the mentioned linearization:

x_ ¼ Ax þ Bu;

xðt0 Þ ¼ x0

y ¼ Cx þ Du; where A, B, C and D are appropriate coefficient matrices of corresponding sizes. The idea of a real-time state estimation also came up in the literature concerning the state observer we apply in the present paper. In [9] real-time estimation process is introduced in purely simulation environment with a numerically stable result. In our present paper, based on one measurement/minute data collection, we realize a short-time follow-up of the state process. Nowadays there is no doubt on the importance to investigate and develop solar thermal heating systems as environmental friendly technology. To our knowledge there is no work in the literature on the design and the application of a state observer for solar thermal heating systems (in short solar heating system in this paper). However, it is an important scope since all state variables are generally not measured in practice (in commercial solar heating systems), though the knowledge of them could be used for energetic evaluation and controlling purposes. In our present work we deal with a mathematical model with four state variables. Two of the state variables are measured and the other two are estimated by a deterministic Luenberger type state observer. Within the discretization naturally determined by the sampling time of the measurements, the state observer corresponds to a continuous time, even time invariant mathematical model. Along with the specified formula of the state observer, the algorithm for its real-time practical application is also introduced as, considering the claimed estimating precision and the discrete time monitoring according to the practical relevance. In the paper the process is also applied to a particular solar heating system producing domestic hot } , Hungary. Simulawater, and monitored by measuring system at the campus of the Szent István University (SZIU), Gödöllo tions and their measurement based checking are used in our investigations. There are several different mathematical models in the literature at different levels of complexity that can be used to investigate the physical process in solar heating systems and to develop them. Duffie and Beckman [10] is a recent comprehensive reference in this topic. Our investigations are based on the mathematical model worked out by Buzás and Farkas [11]. The paper is organized in the following way: Section 2 contains the description of the investigated solar heating system type. In Section 3, the applied mathematical model is introduced. Section 4 is devoted to the construction of the state observer for the model. In Section 5 we give an estimation for the error of the state observer and for the time after which the observer estimates the state with a prescribed precision. In Section 6, the real-time state observer is constructed and the algorithm for the practical use is presented. In this section, based on measurements and the mathematical software package used for the calculations, it is also shown how our state estimation procedure (and the obtained state observer) applies to a particular solar heating system. The measurements and the mathematical model, as well as the estimated state coordinates and their measurements are also compared here. Section 7 contains the obtained results which are discussed in Section 8. Some possible applications are also mentioned in the latter section. In Appendix, from mathematical systems theory, basic concepts and results used in the paper are recalled. 2. Description of the investigated solar heating system type Buzás and Farkas [11] developed a validated linear mathematical model for a quite general kind of solar heating systems, the scheme of which can be seen in Fig. 1. The following time-dependent functions are assigned in the figure: Tc: temperature of the solar collector, °C; Thh: temperature at the hot (collector) side of the heat exchanger, °C; Thc: temperature at the cold (storage) side of the heat exchanger, °C; Ts: temperature of the solar storage, °C; Ic: solar irradiance on collector surface, W/m2; Tca: ambient temperature of the solar collector, °C; Tha: ambient temperature of the heat exchanger, °C; Tsa: ambient temperature of the solar storage, °C; Td: tap water temperature, °C; V_ c : volumetric flow in the collector loop, m3/s; V_ s : volumetric flow in the storage loop, m3/s; V_ l : consumption load of the solar storage, m3/s.

11560

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

Fig. 1. Scheme of the investigated solar heating system type.

Assumptions made here and later: the cooling and delaying effects of the pipes are neglected, the solar collector, both sides of the heat exchanger and the solar storage are completely mixed that is the temperatures (Tc, Thh, Thc, Ts) are homogeneous in these parts of the system. 3. Description of the mathematical model With the notation of the previous section, the linear mathematical model for the solar heating system is the following system of ordinary differential equations:

V_ c Ac g 0 U L Ac T_ c ¼ Ic þ ðT ca  T c Þ þ ðT hh  T c Þ qc c c V c qc cc V c Vc Aa ka qc cc V_ c ekh Ah ðT c  T hh Þ þ c m ðT hc  T hh Þ þ c m 2 ðT ha  T hh Þ Vh Vh h h h h þ qc cc 2 þ qc c c 2 þ qc cc V2h 2 2 2 Aa ka q cs V_ s ekh Ah 2 T_ hc ¼ c m s ðT  T Þ þ ðT  T Þ þ ðT ha  T hc Þ s hc hh hc c h mh c h mh h h þ qs cs V2h þ qs cs V2h þ qs cs V2h 2 2 2

T_ hh ¼ c

h mh

V_ l V_ s As ks T_ s ¼ ðT d  T s Þ þ ðT hc  T s Þ þ ðT  T s Þ; Vs Vs qs cs V s sa where the following constant parameters are used: Ac: collector field surface area, m2; g0: optical efficiency of the solar collector, -; qc: collector fluid density, kg/m3; cc: specific heat capacity of the collector fluid, J/(kgK); Vc: volume of the solar collector field, m3; UL: overall heat loss coefficient of the solar collector, W/(m2K); ch: specific heat capacity of the heat exchanger material, J/(kgK); mh: mass of the empty heat exchanger, kg; Vh: total volume of the heat exchanger, m3; ekh: heat transfer coefficient inside the heat exchanger, W/(m2K) (where e = 1, based on [11]); Ah: surface area of the heat transfer inside the heat exchanger, m2; ka: heat loss coefficient of the heat exchanger to the ambiance, W/(m2K); Aa: surface area of the heat exchanger to the environment, m2; qs: density of water (storage fluid), kg/m3; cs: specific heat capacity of water, J/(kgK); Vs: volume of the solar storage, m3; As: outside surface area of the solar storage, m2; ks: heat loss coefficient of the solar storage to the environment, W/(m2K).

ð1Þ

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

11561

The mathematical model can be handled with the tools of linear systems theory [3], for a more recent reference see e.g. [4]. The state vector in the model consists of four coordinates: Tc, Thh, Thc, Ts, and all of them should be known as functions of time to determine the heat transport and its efficiency at each part of the system. Let us assume that the measurements of one, two or three coordinates of the state vector are available in the solar heating system. This means that only certain state variables are known and from them we want to recover the rest of the state variables as functions of time. (It should be noted that the other time-dependent functions in Section 2 (such as Ic, Tca, etc.) are also measured as inputs at the same time as the measured coordinates of the state vector.) To this end, considering the observability of linear systems, a state observer can be constructed to the considered mathematical model, provided the system is observable, see [5], a more recent reference is [4]. Then the missing state variables can be efficiently estimated by the state observer and can be used, as it is claimed, e.g. to determine the energetical characteristics and efficiencies at different parts of the solar heating system, or to apply some kinds of control on the system. 4. State observer construction for the model of the solar heating system Using the notation and Definition 1 of Appendix, let us rewrite system (1):

2

Tc

2

3

Ic 6 T ca 6 6 u :¼ 6 T ha 6 4 T sa Td

6T 7 6 hh 7 x :¼ 6 7; 4 T hc 5 Ts 2 6 6 6 6 6 A :¼ 6 6 6 6 4 2

_

3 7 7 7 7; 7 5

 qUcLcAVc c  VV cc

V_ c Vc

0

0

qc cc V_ c ch mh V þqc cc 2h 2

qc cc V_ c ekh Ah Aa2ka

ekh Ah c h mh V þqc cc 2h 2

0

0

ekh Ah c h mh V þqs cs 2h 2

0

0

c

Ac g0 q cc V c

6 c 6 6 0 6 B :¼ 6 6 6 0 6 4 0

U L Ac

ch mh V þqc cc 2h 2

A k qs cs V_ s ekh Ah  a2 a ch mh V þqs cs 2h 2

0

0

0

Aa ka 2 ch mh V þqc cc 2h 2

0

0

Aa ka 2 ch mh V þqc cc 2h 2

0

0

0

qc cc V c

As ks

qs cs V s

0

_

qs cs V s ch mh V þqs cs 2h 2 V_ l V_ s Vs

V_ s Vs

3

 qAcs skVs s

7 7 7 7 7 7; 7 7 7 5

s

3

7 7 07 7 7; 7 07 7 5

V_ l Vs

where x is considered the state of the system (its coordinates are the state variables) the other quantities of Fig. 1 (Ic, Tca, etc.) will be considered constant on a given time interval (see Section 6). Let us assume that V_ c ; V_ s and V_ l are constant in time so A and B are also constant. Hence the observation system

x_ ¼ Ax þ Bu; y ¼ Cx;

ð2Þ

is considered, where observation matrix C may be any projection matrix, corresponding to the observation of one or several state coordinates, e.g.

½ 0 1 0 0 ;



1 0 0

0

0 0 1 0



2

1 0 0

0

3

6 7 ; 4 0 0 1 0 5: 0 0 0 1

Let us seek a state observer of the form as in Definition 2 of Appendix:

z_ ¼ Fz þ Bu þ Hy:

ð3Þ

Thus we get the so-called error dynamics:

_ ¼ Fz þ Bu þ HCx  Ax  Bu: e_ ¼ ðz_  xÞ

ð4Þ

Let the following relation hold between F and H: F = A  HC. Then from (4) we get

e_ ¼ ðA  HCÞe:

ð5Þ

11562

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

From Definition 2 of Appendix, we can see that system (3) is a state observer for (2) if system (5) is asymptotically stable, that is, the real parts of all eigenvalues of A  HC are negative. For the determination of such H, under certain conditions, we can even prescribe pair-wise conjugate eigenvalues k1 ; . . . ; kn of A  HC with negative real parts (see Remark 1 of Appendix). In other terms, we want to solve a pole assignment problem for the pair (A, C) (see Definition 3 and Remark 1 of Appendix.) According to Theorem 1 of Appendix, this can be done if the so-called Kalman matrix ½C 0 ; A0 C 0 ; : . . . :; ðA0 Þn1 C 0  has full rank (where prime denotes the transpose of the matrix). The latter condition is known to be equivalent to observability of system (A, B, C) which means that the state process x can always be recovered uniquely from the observation y. 5. Estimation error and minimal threshold time for an error bound Let us consider the error dynamics (5). Here and later we consider the only case when the eigen values of A  HC are given different, strictly negative real numbers. Then the solution of system (5) is of the form

eðtÞ ¼

n X ci eki t si

ðt 2 Rþ0 Þ;

ð6Þ

i¼1

2

3 si1 6 : 7 6 7 7 where ki ði ¼ 1; . . . ; nÞ are the eigenvalues of A  HC, si :¼ 6 6 : 7ði ¼ 1; . . . ; nÞ are the corresponding eigenvectors, with 4 : 5 sin |s1| =    = |sn| = 1, and ci 2 R (i = 1, . . ., n) depend on e(0). In the following let us consider such error vectors e(0), the absolute values of which are less than or equal to the maximal possible absolute error, R, given by the current physical arrangements. (See Section 6 for the study of a particular application.) Let us take an error ball with radius R. (See Fig. 2.) Let us fix an initial error vector e1(0) in the space Rn, starting from the origin, with maximal length, R and perpendicular to the hyperplane S1 (generated by eigenvectors s2, ..., sn). (Two vectors fulfil these conditions.) Let c1 denote the coefficient of eigenvector s1 in formula (6). By fixing e2(0),...,en(0), with similar constructions c2 ; . . . ; cn can be also determined. Now, for any arbitrary vector e(0) in the error ball (with length less than or equal to R) for all coefficients of e(0) according to the form (6) we have jci j 6 jci j (i = 1, ..., n). Finally the coordinates of the error vector, e(t), in the canonical orthonormed basis of Rn, can be estimated from above by the following term strictly monotonically decreasing with time:

jek ðtÞj 6 jeðtÞj < jc1 s1 ek1 t j þ    þ jcn sn ekn t j < jc1 jjs1 jek1 t þ    þ jcn jjsn jekn t ¼ jc1 jek1 t þ    þ jcn jekn t ðk ¼ 1; . . . ; nÞ:

ð7Þ

Based on estimation (7) we can calculate the time (s) from which each absolute error, |ek(t)| (k = 1, ..., n) is surely less than a prescribed error bound (e). The latter holds, whenever the right hand side of (7) is less than e. Now let us determine the state estimate vector z. The coordinates of z estimate precisely enough the coordinates of x from time s on. The calculation of the coordinates of z, that is the estimation of the coordinates of x, is performed based on real-time measurements, namely by the monitoring of vector u, V_ c ; V_ s ; V_ l and the measured coordinates of vector x .

Fig. 2. Initial error ball.

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

11563

The initial values zi(0) of the coordinates of z, corresponding to the measured coordinates of x are the measured coordinates of x(0) themselves. For the initial values of the other coordinates of z(0), arbitrary, physically relevant values can be taken so that the absolute value of the initial error |e(0)| is not greater than R of Fig. 2. Since within the prescribed or given time interval of length Dtj = tj  tj1 (j = 1, . . ., T), between the (j  1)st and the jth sampling time of the measurements (where t0 = 0, t0 < t1 <    < tT), the monitored variables are taken constant, the method discussed above can be applied within every time interval, and a real-time estimation of the missing state variables can be carried out by z with the prescribed precision (e) as detailed in the next section. 6. Construction of real-time state observer 6.1. Description of the algorithm For our four-dimensional model (n = 4), starting from t0 = 0, at every measured time point tj (j = 1, ..., T), z(tj) can be determined by z(tj1) according to the following cases. Fig. 3 sketches the algorithm of the real-time state observer design according to the method below. j k 1. If rank C 0 ; A0 C 0 ; . . . ; ðA0 Þn1 C 0 ¼ 4 at tj1 (i.e. the Kalman rank condition holds), then for system (A, B, C), a state observer can be constructed, according to the measurements at tj1, as written in Section 4, so that the eigenvalues of A  HC are arbitrary different negative real numbers prescribed by us. Then the value of s is calculated to the prescribed error bound (e) according to Section 5. If s < Dtj, the adequate state observer can be obtained from (3).

Fig. 3. Algorithm of real-time state observer design.

11564

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

If s P Dt j , then to speed up the convergence, each eigenvalue must be increased in module. In our algorithm, we multiply each prescribed eigenvalue by 1.1 (several times, if necessary). The resulting eigenvalues are denoted by ki , see Fig. 3. This process should be repeated till s < Dtj is fulfilled. Then with initial value z(tj1), z(tj) can be determined from equation (3). 2. If rank½C 0 ; A0 C 0 ; . . . ; ðA0 Þn1 C 0  < 4 at tj1, (i.e. the Kalman rank condition does not hold), then a state observer is not constructed, we just set z(tj) :¼ z(tj1). The process stops when all z(tj)  s (j = 1, . . ., T) are determined. For a graphic representation we will fit a polygon to these points where the state estimates are the ‘‘best’’. u; V_ c ; V_ s ; V_ l and the measured coordinates of x are taken as piecewise constant functions according to the time steps of the measurements. To each observation C there corresponds an algorithm of the above type. In the special case of observing only Tc, we have   1 0 0 0 when both are observed. Considering the definition C ¼ ½ 1 0 0 0 , for Ts only C ¼ ½ 0 0 0 1 , and C ¼ 0 0 0 1 of matrix A, (in Section 4), it can be easily seen that if the pumps are off, the Kalman rank condition does not hold for any observation above. On the contrary in our investigations, if the pumps were on, we always experienced that the rank condition was fulfilled for any of the observations, thus the state observer could be always constructed. In accordance with the customary operation of the solar heating system, both pumps are on or off simultaneously every time. The consideration of the above observations is justified by the fact that the collector outlet temperature, Tc, and the storage temperature, Ts, are available in commercial solar heating systems since these values are needed for the controller unit. 6.2. Observer design for a commercial solar heating system In this section the above methodology is applied to a solar domestic water heating system with flat plate collector field and external heat exchanger, see [12]. For a given date we have got the measurements of all state variables (Tc, Thc, Thh, Ts) for the whole day. Then with the observed state variables (Tc and Ts which are considered measured now), we calculate the state estimate vector z(t) for the same day. After that, all z(t) will be compared with the measurements. All calculations were done in Matlab 7.1 simulation environment [13]. First of all we validate the mathematical model (1) for our particular solar heating system by comparing measured state variables with calculated ones by the model for a given date. Parameter values of the solar heating system: g0 = 0.74, Vc = 0.027 m3, mh = 37 kg, Vh = 0.005 m3, ekh = 2461.5 W/(m2K) and Ah = 2 m2 are from the technical documentation of the system. UL = 7 W/(m2K) is from recommendation in the literature [14]. ka = 5 W/(m2K) and ks = 1 W/(m2K) are estimated values based on measured data. Ac = 33.3 m2; qc = 1034 kg/m3; cc = 3623 J/(kgK); ch = 464.76 J/(kgK); Aa = 0.24 m2; qs = 1000 kg/m3; cs = 4200 J/(kgK); Vs = 2 m3; As = 4 m2. The pumps operate simultaneously and their flow rates, V_ c and V_ s are measured. V_ l and all coordinates of both vectors u and x are measured as a function of time, too. More precisely, for technical reason, the outlet pipe temperatures from the solar collector field and from both sides of the heat exchanger are measured and taken as measured values of Tc, Thh, and Thc. The solar storage temperature Ts is measured inside the storage. The monitored values of Thh and Thc serve only for comparison purposes. The measurements are carried out once in every minute. The measured day is 3rd April, 2011.

Fig. 4. x1 and Tc,meas on 3rd April, 2011.

11565

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

Fig. 5. x2 and Thh,meas on 3rd April, 2011.

Fig. 6. x3 and Thc,meas on 3rd April, 2011.

Fig. 7. x4 and Ts,meas on 3rd April, 2011.

2

3 T c;meas ð0Þ 6 T hh;meas ð0Þ 7 7 The state vector x of model (2) is also calculated for the same day, with the measured initial value xð0Þ ¼ 6 4 T hc;meas ð0Þ 5. T s;meas ð0Þ Figs. 4–7 compares the measured and modelled state variables. The operating state of the pumps (on/off) can be also seen in each figure. Fig. 4 shows the modelled temperature of the solar collector field, x1 (with dashed line) and the measured outlet temperature of the solar collector field, Tc,meas.

11566

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

Fig. 5 shows the modelled temperature of the hot side of the heat exchanger, x2 (with dashed line) and the measured outlet temperature of the hot side of the heat exchanger, Thh,meas. Fig. 6 shows the modelled temperature of the cold side of the heat exchanger, x3 (with dashed line) and the measured outlet temperature of the cold side of the heat exchanger, Thc,meas. Fig. 7 shows the modelled temperature of the solar storage, x4 (with dashed line) and the measured temperature of the solar storage, Ts,meas. The average of the absolute difference between the modelled and measured temperatures for the given date is 3.7 °C in case of the solar collector field, 2.4 °C in case of the hot side of the heat exchanger, 2.9 °C in case of the cold side of the heat exchanger and 0.3 °C in case of the solar storage. In proportion to the difference between the maximal and minimal measured temperature values, the time average of the absolute error is 6.6% in case of the collector field, 6.2% in case of the hot side of the heat exchanger, 7.9% in case of the cold side of the heat exchanger and 3.9% in case of the solar storage. Based on Fig. 4–7 and the above evaluation, we state that the model tracks the physical processes characteristically right with an acceptable accuracy in view of several engineering aims. So we accept the used mathematical model (1), noting that on a certain time interval the largest discrepancy between the measurements and calculations occurs in case of the solar collector field. The model cannot take into account the random effect of sky cloudiness which is likely to cause this error. Let the accuracy prescription be e = 0.1 °C, as defined in Section 5, that is

jek ðtÞj 6 0:1  C

ðk ¼ 1; . . . ; 4Þ:

ð8Þ

The state observer is chosen such that from the state vector the collector and the solar storage temperatures are observed     Tc 1 0 0 0 (measured). y ¼ . With the notation of observation system (2), we have C ¼ , providing Ts 0 0 0 1 0 0 0 2 0 0 3 0 0 _ _ rank½C ; A C ; ðA Þ C ; ðA Þ C  ¼ 4, if V c – 0 and V s – 0 . Therefore, a state observer can be constructed in this case. In our estimation process, in addition to the monitored state variables Tc and Ts, from the state observer we get estimations for the unmeasured variables, Thh and Thc at the end of each time interval Dtj . The radius R of the error ball (see Section 5) now by physical considerations is

2 3  0    6 100 7 6 7 R ¼ 6 7; 4 100 5    0  where the first and the last coordinates are 0 because Tc and Ts are observed without estimation error. The other two coordinates are the differences between the highest and lowest thinkable values of the corresponding temperatures. These values are the maximum absolute error that we may make at any estimation. Let the prescribed initial value of z be

2

T c ð0Þ

3

6 T ð0Þ 7 7 6 c zð0Þ :¼ 6 7; 4 T s ð0Þ 5 T s ð0Þ where the first and the last coordinates are directly measured data which are also used for the initial estimation of the second and third coordinates. In our calculations the prefixed eigenvalues are k1 :¼ 1; k2 :¼ 0:9; k3 :¼ 0:8; k4 :¼ 0:7 which turned out to be suitable for each time interval, since we always obtained s < Dtj = 1 min (j = 1, ..., T). The matrix A  HC is determined at every time step where the rank condition holds and the following state observer is obtained.

3 Ic 7 6 6 T ca 7   7 6 Tc 7 þ H z_ ¼ ðA  HCÞz þ Bu þ HCx ¼ ðA  HCÞz þ B6 : T 6 ha 7 Ts 7 6 4 T sa 5 Td 2

7. Results

ð9Þ

2

3 2 3 T c;meas z1 6 T hh;meas 7 6 z2 7 7 6 7 By the measuring and the estimation process, the measured and the estimated state vectors, 6 4 T hc;meas 5 and z ¼ 4 z3 5, have been specified at every measurement instant for the particular date, 3rd April, 2011. T s;meas z4

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

11567

Fig. 8. z2, Thh,meas and the operating state of the pumps on 3rd April, 2011, Tc and Ts are observed.

Fig. 9. z3, Thc,meas and the operating state of the pumps on 3rd April, 2011, Tc and Ts are observed.

Figs. 8 and 9 compare the measured and the estimated state variables. The operating state of the pumps (on/off) can be also seen in each figure. (The estimated values of Tc and Ts are not considered since they are directly measured.) Fig. 8 shows the estimated temperature of the hot side of the heat exchanger z2 and the measured outlet temperature of the hot side of the heat exchanger Thh,meas. Fig. 9 shows the estimated temperature of the cold side of the heat exchanger z3 and the measured outlet temperature of the cold side of the heat exchanger Thc,meas. It can be seen from the figures that the estimation accuracy is more acceptable when the pumps are switched on. This is not surprising since the state observer can be constructed in this case. Hence we make the error evaluation for this time, too. The average of the absolute difference between the estimated and measured temperatures for the time when the pumps are switched on is 2.9 °C in case of the hot side of the heat exchanger and 1.8 °C in case of the cold side of the heat exchanger. In proportion to the difference between the maximal and minimal measured temperature values for the given date, the average of the absolute error for the time when the pumps are switched on is 7.5% in case of the hot side of the heat exchanger and 5.2% in case of the cold side of the heat exchanger. 8. Discussion and outlook Since the applied mathematical model (1) is based on physical foundations (white-box model) it can be required to approximate well the real operation of the investigated solar heating system. This expectation is fulfilled as it can be stated that, considering the comparison with measurements in Section 6.2 and taking the aims of the model (studying the main processes and developing the efficiency of solar heating systems), it estimates the real processes characteristically right and, for the general engineering practice with convenient accuracy. Namely for studying the main physical processes to a first simple (linear) approximation and, based on the latter, developing solar heating systems.

11568

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

The apparent difference between the measured and the modelled graphs is mostly because of the fact that the measurements have been made not inside but on the outlet pipes of the collector field and both sides of the heat exchanger. Accordingly, it should be noted that in case of switched-on pumps, the fit of the model to the measurements is better. The used model is linear and deterministic. A richer model can be obtained by building in random processes (stochastic description of sky cloudiness). Based on the calculations and measurements, it can be stated that the state estimate vector z estimates the real data more reasonably if the pumps are on for a while permanently. This is not surprising, since the existence of a state observer is guaranteed by the rank condition only in case of switched on pumps, and as mentioned above, the model fits better, implying a better performance of the observer. (This experience is in accordance with the intuitive expectation that in case of switched off pumps, the state of the system cannot be uniquely reconstructed from the measured state variables. E.g. in case of no flow in the system, Tc and Ts can be modulated by solar irradiation and consumption load independently from each other and from Thh and Thc .) If the pumps are off, the coordinates of z are constant, so it should be noted that the use of our state observer is recommended only in the time periods of switched on pumps, e.g. for the applications below. Although the state observer has its limits in precision, its simplicity (linearity) is an important advantage and it enables us to develop a real-time estimate of the unmeasured state variables of solar heating systems that can be used for several purposes in engineering practice. We made measurements once in every minute. Increased measurement frequency should cause better state observer performance in accuracy. As for a possible nonlinear model, we think that a similar estimation method can be developed near an equilibrium. As for the possible applications of our methodology, in commercially available systems Tc and Ts are generally measured but other variables are not. The estimation of Thh and Thc by the state observer makes the energetical and efficiency calculations of different parts of a solar heating system possible. Another possible application concerns the control of solar heating systems. The generally used ordinary control method switches on the pumps if the measured temperature difference of Tc  Ts is higher than a prefixed constant value DTon. This value is set relatively high in order to surely avoid the cooling down of the solar storage by the collector field which results in some omission of the available solar energy (irradiation.) By measuring Tc and Ts, Thc can be real-time estimated with some assumed values of V_ c and V_ s by the here developed state observer. This estimated temperature can be taken as inlet temperature to the solar storage. If this value is high enough to efficiently switch on the pumps, considering the amount and the cost of the used electric energy and saved gas energy, then the pumps are switched on. Alternatively, with the help of the state observer and such a control method, in case of changing flow rate mode of the pumps, V_ c and V_ s could be set at every instant so that it is worth to switch them on even if the saving is minimal. A higher rate of the available solar irradiation may be utilized if this energetically-based control operates instead of the generally used ordinary control, or if the energetically-based and the ordinary controls are connected with the logical OR operator. The latter means that the pumps are on if either control orders the pumps to switch on. For the detailed analysis of a control method that works with the prediction of the heat exchanger outlet temperature, see [15]. Acknowledgements The authors say special thanks to Dr. János Buzás, Prof. István Farkas and the Physics and Process Control Department in the SZIU for the contribution to this work and for giving the possibility to carry out measurements on the considered solar heating system. Appendix For the reader’s convenience, from the standard mathematical systems theory ([3], and a more recent reference is [4]), we recall some basic concepts and theorems used in the paper. Definition 1. Given n, m, p 2 N, A 2 Rnxn, B 2 Rnxm, C 2 Rpxn, the triple (A, B, C) is called a (linear time-invariant) controlobservation system with the following interpretation: For fixed initial state x(0) and piecewise continuous input (control) function u, the state process x of the system is the corresponding solution of

x_ ¼ Ax þ Bu

ðA1Þ

and instead of the state x its transformed

y ¼ Cx;

ðA2Þ 1

p

is observed, where y:R ? R is called observation.In the paper constant inputs are considered. Definition 2. Let (A, B, C) be a control-observation system according to Definition 1.With F e Rnxn, H e Rpxn,

z_ ¼ Fz þ Bu þ Hy;

ðA3Þ

R. Kicsiny, Z. Varga / Applied Mathematics and Computation 218 (2012) 11558–11569

11569

is called an observer system (state observer) for system (A, B, C) if for any x(0), z(0) initial conditions and u e U, for the solutions of (A1) and (A3), limt!1 ðzðtÞ  xðtÞÞ ¼ 0 holds.Z: R1 ? Rn is called state estimate, and e := z  x is the estimation error. Definition 3. Let M1 e Rnn, M2 e Rpn, ai e R, i = 0, . . ., n  1. If for an arbitrary polynomial pðkÞ ¼ kn þ an1 kn1 þ    þ a0 , there exists D e Rnp so that detðkI  ðM 1 þ DM2 ÞÞ ¼ pðkÞ, where I e Rnn is the identity matrix, then we also say that the pole assignment problem for the pair (M1, M2) can be solved. Remark 1. Since pðkÞ is arbitrary, Definition 3 means that for any pair-wise conjugate complex numbers k1 ; . . . ; kn there exists D e Rnxp such that k1 ; . . . ; kn are the eigenvalues of M1 + DM2. Theorem 1. With the notation of Definition 3, rank½M 02 ; M 01 M 02 ; : . . . ; ðM 01 Þn1 M 02  ¼ n implies that the pole assignment problem for the pair (M1, M2) can be solved. References [1] N.M. Abdel-Jabbar, R.Y. Jumah, M.Q. Al-Haj Ali, State estimation and state feedback control for continuous fluidized bed dryers, Journal of Food Engineering 70 (2005) 197–203. [2] M. Perrier, S. Feyo de Azevedo, E.C. Ferreira, D. Dochain, Tuning of observer-based estimators: theory and application to the on-line estimation of kinetic parameters, Control Engineering Practice 8 (2000) 377–388. [3] R. Kalman, M. Arbib, P. Falb, Topics in Mathematical System Theory, McGraw-Hill, New York, 1969. [4] S.H. Zak, Systems and Control, Oxford University Press, New York-Oxford, 2003. [5] D.G. Luenberger, Observing the state of a linear system, IEEE Transactions on Military Electronics 8 (2) (1964) 74–80. [6] F. van der Heijden, R.P.W. Duin, D. de Ridder, D.M.J. Tax, Classification, parameter estimation and state estimation, An Engineering Approach using MATLABÒ, John Wiley & Sons, Ltd, 2004. [7] J. O’Reilly, Observers for Linear Systems, Mathematics in Science and Engineering, Academic Press Inc., London, 1983. [8] C. Nguyen, T. Lee, Design of a state estimator for a class of time-varying multivariable systems, IEEE Transactions on Automatic Control 30 (2) (1985) 179–182. [9] T. Grootscholten, K. Keesman, P. Lens, Modeling and on-line estimation of zinc sulphide precipitation in a continuously stirred tank reactor, Separation and Purification Technology 63 (2008) 654–660. [10] J.A. Duffie, W.A. Beckman, Solar Engineering of Thermal Processes, John Wiley and Sons, New York, 2006. [11] J. Buzás, I. Farkas, Solar domestic hot water system simulation using block-oriented software, in: The 3rd ISES-Europe Solar World Congress (Eurosun 2000), Copenhagen, Denmark, CD-ROM Proceedings, (2000) p. 9. [12] I. Farkas, J. Buzás, A. Lágymányosi, I. Kalmár, E. Kaboldy, L. Nagy, A combined solar hot water system for the use of swimming pool and kindergarten operation, Energy and the environment, in: B. Frankovic(Ed.), Croatian Solar Energy Association, Opatija, vol. I. (2000) pp. 81–88. [13] D.M. Etter, D. Kuncicky, H. Moore, Introduction to MATLAB 7, Springer, 2004. [14] B. Bourges, European Simplified Methods For Active Solar System Design, Kluwer Academic Publishers for CEC, Available from the Energy Research Group, University Collage, Dublin, 1991. [15] R. Kicsiny, Development of an energetically-based control for solar thermal heating systems, Review of Faculty of Engineering, Analecta Technica Szegediensia, Szeged, (2009) pp. 50–57.