# Central Heating System Constrained Control with Input Delay Based on Neural Networks

Authors: Hongwei Wang; Fangwen Tu; Guohui Feng; Xin Ao
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

An output constrained control with input delay is proposed for a central heating system. Due to the delay of signal transmission and valves opening time, an input delay is considered into the system and an auxiliary system is employed to handle this issue by converting the delayed input into a delay-free one. Moreover, to ensure the output supply water temperature within a limited range, Barrier Lyapunov algorithm is involved to achieve desired control accuracy. Finally, external disturbance and model uncertainty are incorporated into the dynamic system and neural networks (NN) are trained in an online fashion for the compensation. The stability of the control system is guaranteed through rigorous Lyapunov analysis and the excellent control performance over traditional PID control is demonstrated via numerical simulation study.


---

## Body

## 1. Introduction

Currently, central heating system is playing an essential role to human’s daily life and the boast of economy in everywhere all over the world. Due to the development of technology, advanced and intelligent central heating system has been widely researched and deployed to the municipal service system. This kind of system has the advantages of energy saving and contributes a lot to the less-pollution environment. Typically, this system measures the indoor/outdoor temperature and other information to control the valve opening degree in primary pip network which indirectly controls the supply water temperature in secondary network and provides heating service to the users.In literature, many researches have been done to estimate or calculate the heating load which subsequently can help to adjust the heating control scheme [1–3]. However, it is rare research which focuses on the automatic control of central heating system. The representative works that can be tracked include [4] in which an indirect connection of district heating system architecture is built to achieve an optimal control. In [5], a regional automatic control strategy is proposed to designate temperature monitor and valve on the user end for the heating quality enhancement. A load prediction based multiple constrained control is put forward in [6] to handle the accuracy of water temperature supply and input water volume constraints. However, none of the above mentioned research considers the input delay due to the signal transmission inside the system and the valve operating time. It brings additional challenges for the control design since the time-lag will easily cause the instability of the whole system. For linear system, many researches in the literature can be found to cope with this problem [7, 8]. But for the complicated nonlinear system, it becomes more difficult to design a proper control. In [9], a virtual observer is developed to act as an auxiliary system that can convert the input delay system into a nondelay one. The target control system has a multi-input and multioutput (MIMO) nonlinear structure. A predictor of an uncertain nonlinear system is developed for the compensation of delay in the input to upgrade the control performance of a PID control. In this paper, inspired by the aforementioned researches, we will convert the original input delay heating system into a delay-free one similar to the Artstein model [10] together with the output constrained control which will be introduced in the following contents. This design will enhance to the stability of the delay system via providing designed accumulated compensation. The excellent performance will be verified in the simulation section.Apart from input delay, the accuracy of supply water temperature in the secondary network is also of significance to the heating quality. Therefore, in this work, the output tracking error constraints are novelly taken into consideration for the control design. Some of the existing researches offer good results such as the ones using artificial potential field [11], prescribed performance control [12, 13], model predictive control [14] as well as reference governor [15]. Additionally, Barrier Lyapunov Function (BLF) is introduced in [16–18]. This method requires less restrictive initial conditions and does not rely on explicit system solution. Due to its advantage, BLF is tentatively incorporated into the target heating system with some customized modifications to guarantee the limited tracking error of the supply water temperature. Finally, the external disturbance, model uncertainty error, and other unknown factors are modeled as a disturbance term complemented into the system dynamic model. Different from a robust control [11, 19], in this paper, we design an NN [20] observer to estimate the disturbance to achieve more accurate control and attenuate the noise. The overall system diagram of the proposed central heating control system is summarized in Figure 1.Figure 1
Control system architecture diagram.The contributions of this paper are threefold:(i)
Different from existing heating system control research, input delay due to the action of valve and signal transmission is novelly considered in the nonlinear system. The delay system is converted into a delay free one through an axillary updating system for the compensation.(ii)
The input delay control is proposed in combination with output constraints that can regulate the tracking accuracy of the supply water temperature in the secondary loop with specifically designed BLF method.(iii)
The external disturbance of the system and the model error are considered and handled with an NN approximator. Weights updating law is developed with rigorous mathematical verification. And the outstanding performance is demonstrated through simulation study.The organization of this paper is as follows. Section2 introduces the birdview of the overall system and dynamic model. The detailed mathematical deduction of the proposed control is presented in Section 3. Section 4 conducts a numerical simulation study to demonstrate the efficiency of the proposed scheme. Finally, some concluding remarks are reached in Section 5.


## 2. System Overview and Dynamic Model

In this work, the target heat supply system is an indirectly connected central heating system with excellent efficiency and energy conservation capability. The whole picture of this system can be depicted as Figure2 where primary and secondary network are connected together via heat exchanger for the purpose of heat transformation. The heat transformation process follows the heat balance principle. In the following contents, we focus on the secondary network since this loop influences directly the users’ well-being and heating quality. The mathematical model is put forward as [21](1)Q1=Q2=Q3(2)Q1=qvVwTn-Tw(3)Q2=a·FwTg2+Th22-Tn1+b(4)Q3=G2cwTg2-Th2where Q1 represents design building space-heating load, Q2 denotes heat from radiator, and Q3 is the heat provided to the users through heat supply loop. qv, VW, a, Fw, G2, and cw stand for the heating index, building volume, scalar parameter, heat radiating area of the radiator, design flow of secondary network, and water specific heat capacity, respectively. b is fixed as 0 for hot water floor radiation heating. Tn and Tw denote indoor and outdoor temperature. Tg2 and Th2 represent supply and return water temperature in secondary network. Following the heat balance principle, the mentioned three heat loads Q1, Q2, and Q3 are supposed to be equal to each other as presented in (1).Figure 2
System architecture of indirectly connected district heating system.Assumption 1.
Secondary network usually has 3-5 % heat loss. For simplicity purpose, in this work, the loss is neglected.Assumption 2.
In this paper, the daily real-time heating load and desired secondary network supply water temperatureTg2d are known.


## 3. NN Based Adaptive Control with Input and Output Constraints

### 3.1. Preliminaries on NN Approximation

Radial Basis Function Neural Network (RBFNN) has been demonstrated to have outstanding function approximation and learning capabilities. If we have a continuous nonlinear functionϕ(z):Rn→Rm defined on a compact set Ωz, it can be approximated by RBFNN with following expression:(5)ϕNNZ=WTSNNZwhere Z∈Ωz⊂Rn represents the input vector. W∈Rp×m stands for the adjustable weight matrix and p is the number of neuron. SNN(Z)=[s1(Z),…,sp(Z)]T is the basis function vector, where si(Z)=exp⁡[-(Z-μi)T(Z-μi)/ψi2],i=1,2,…,p. μi=[μi1,μi2,…,μin] are the centers of respective field and ψi denotes Gaussian function width. Concretely, (5) can be described as(6)ϕZ=W∗TSNNZ+εZwhere W∗ is the ideal weight and ε(Z)∈Rm denotes the approximation error and is subject to ε(Z)≤σ, where σ should be a positive constant. The ideal weight matrix W∗ can be defined as(7)W∗≔arg⁡minW∈Rp×m⁡supZ∈ΩZ⁡ϕZ-WSNNZIt should be noted that W∗ is only an “artificial" quantity for the purpose of analysis. In control design, the ideal values will be estimated by W with developed updating law [20]. Moreover, NN approximation error ε(Z) indicates the minimum deviation between optimal approximation solution, i.e., W∗SNN(Z) and the unknown function ϕ(Z). It has been proved that approximation with neural networks can theoretically have any expected accuracy if the neuron number is large enough; i.e., ε(Z) is able to be arbitrarily small if p is sufficiently large [20].


### 3.2. Adaptive Control with Multiple Constraints

#### 3.2.1. Mathematical Model Formulation

If we are already given the supply water temperatureTg2 in secondary network, our main aim is to control the valve opening degree in primary network to guarantee the supply water temperature following the desired values. Concretely, the given temperature value Tg2 acts as the desired signal symbol as Tg2d. The transfer function in this heat exchange process can be expressed as [22](8)G1s=KpT1s+1T2s+1e-τ1swhere Kp denotes a scalar coefficient, T1 and T2 represent the constant inertial time, and τ1 is the lagging time. This mathematical model depicts the transformation from system input Vg1, i.e., opening degree of the valve in primary pipe network, to the system output Tg2, i.e., actual water temperature in secondary pipe network. For the purpose of control design, we convert the model above into a time-domain formulation.(9)T1T2KpT¨g2+T1+T2KpT˙g2+Tg2Kp=Vg1t-τ1Equation (9) can be concisely written as(10)MT¨g2+CT˙g2+gTg2=Vg1t-τ1where M=T1T2/Kp, C=T1+T2/Kp, and g(Tg2)=Tg2/Kp. In order to extend the generalization of our control, we consider the following remarks and assumptions.Remark 3.
Even though the parametersKp, T1, and T2 in [22] are assumed as constants, in practical application, it is apparent that the heat exchange process is time-varying. In other words, the aforementioned coefficients should be state-dependent and time-varying. In this regard, M, C, and g(Tg2) can be more accurately considered as M(T˙g2,Tg2,t), C(T˙g2,Tg2,t), and g(Tg2,t).Remark 4.
In actual system, external disturbance as well as the model uncertainty should be considered in model (10). For simplicity purpose, they are merged into single term as d(Tdg2,Tg2,t).Besides the remarks and assumption made above, one more variableTdg2=T˙g2 is considered which is actually denotes the temperature variation in secondary loop. Thus, the system model of (10) can be modified as(11)T˙g2=Tdg2MTdg2,Tg2,tT˙dg2+CTdg2,Tg2,tTdg2+gTg2,t=Vg1t-τ1+dTdg2,Tg2,tAssumption 5.
In many researches, only single primary/secondary pipe network is considered; i.e., aforementioned model is defined as an SISO system. In this paper, system is extended to be an MIMO system which can produce multiple water temperatures in secondary pipe network to satisfy different operating requirements.Assumption 6.
The inertia matrixM is invertible and M-1 is bounded. The upper bound will be denoted as M-1∞≤M-1¯, where M-1¯∈R is a positive constant bound.With Assumption5, the input and output variables become a second order vector Tg2=[Tg21,Tg22]T and Vg1=[Vg11,Vg12]T. Additionally, since the disturbance term d(Tdg2,Tg2,t) is going to be estimated with an NN observer, it is neglected temporarily. The resultant system model can be concisely expressed as follows. Herein, we slightly abuse the use of M(Tdg2,Tg2,t), C(Tdg2,Tg2,t) to M, C(12)T˙g2=Tdg2MT˙dg2+CTdg2+g=Vg1t-τ1


#### 3.2.2. Adaptive Control Design without Disturbance

In this paper, the output constraints, i.e., the tracked secondary loop water temperature errors, will also be considered to regulate within a desired range. In order to handle the output constraint problem, Symmetry Barrier Lyapunov Function (SBLF) [23] is employed together with a backstepping approach for the control development.Step 1.
Denote(13)z1=Tg2d-Tg2and(14)z2=αc-Tdg2where the desired trajectory satisfies Tg2d,T˙g2d∈L∞. αc refers to stabilizing function. Select a positive definite and C1 continuous SBLF candidate as(15)V1=12logNbTIxNbNbTIxNb-z1TIxz1+12logNbTIyNbNbTIyNb-z1TIyz1where(16)Ix=1000Iy=0001Nb∈R2×1 denotes the tracking error constraint such that z1≤Nb should be satisfied through the control.Remark 7.
In practical application, the initial conditions of water temperature and temperature variation are the same with the desired values. Hence,z10<Nb should be satisfied.Time derivative ofV1 gives(17)V˙1=z1TIxz˙1NbTIxNb-z1TIxz1+z1TIyz˙1NbTIyNb-z1TIyz1Differentiating z1 with respect to time gives(18)z˙1=T˙g2d-αc-z2Substituting (18) into (17), we obtain(19)V˙1=z1TIxT˙g2d-αc-z2NbTIxNb-z1TIxz1+z1TIyT˙g2d-αc-z2NbTIyNb-z1TIyz1Design the stabling function αc as(20)αc=T˙g2d+NbTNb-z1Tz1K1z1Substituteing (20) into (19), the following can be derived:(21)V˙1=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1Step 2.
Invoke an auxiliary stateS∈R2×1 to cope with the input delay which has the expression following formulation:(22)S=z2-M-1∫t-τ1tVg1θdθ-zfwhere zf∈R2×1 is developed with the following adaptive law.(23)z˙f=K2S-Γ1z2-Θzfwhere K2,Γ1,Θ∈R2×2 are positive coefficients. Multiplying both sides of (22) with M and involving Ms=C(Tdg2,Tg2,t)Tdg2+g(Tg2,t), the derivative of MS gives(24)MS˙=Mz˙2-Vg1t+Vg1t-τ1-z˙f=Mα˙c+CTdg2,Tg2,tTdg2+gTg2,t-Vg1t-K2S+Θzf+Γ1z2=Mα˙c+Ms+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2where without loss of generality M is considered as a constant matrix since the inertial time constants change very slowly. The time-varying part can be added to the disturbance term d. Moreover, Nc is designed as follows with the consideration of Mean Value Theorem [24]:(25)Nc=Θzf+Γ1z2+K2z2+ST+S˙Tz2Nc≤N¯czszswhere the bounding function N¯c(zs) denotes a globally positive function. zs is defined as zs=[z1T,z2T,ST,zτT,zfT]T, where zτ∈R2×1 refers to(26)zτ=Vg1t-Vg1t-τ1=∫t-τ1tV˙g1θdθAfter the introduction of auxiliary state S, the input-delay system is transferred into a delay-free system as expressed in (24). In order to prove the stability of the resultant system, a quadratic form Lyapunov-Krasovskii candidate function can be designated as [25](27)V2=V1+12z2Tz2+12STMS+12zfTzf+υ∫t-τ1t∫wtV˙g1θ2dθdw
DifferentiateV2 and incorporate (21), (22), (23), and (24); this yields(28)V˙2=V˙1+z2Tz˙2+STMS˙+zfz˙f+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1+z2TS˙-M-1Vg1t-τ1-Vg1t+K2S-Θzf-Γ1z2+STMα˙c+Ms+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2+zfTK2S-zfTΘzf-zfTΓ1z2+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1-z2TΓ1z2+z2TM-1zτ-z2TΓ1+Izf-STK2S+STMα˙c+Ms+Nc-Vg1θ+zfTK2S-zfTΘzf+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθThe control law can be developed as follows:(29)Vg1t=Mα˙c+Ms+K2zf+ST+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1+NbTIxNbz1TK1z1NbTIxNb-z1TIxz1+NbTIyNbz1TK1z1NbTIyNb-z1TIyz1Substituting (29) into (28) and invoking (25) and Assumption 6, we have(30)V˙2=-2z1TK1z1-z2TΓ1z2-STK2S-zfTΘzf+z2TM-1zτ-z2TΓ1+Izf+STNc+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1-NbTIyNbz1TK1z1NbTIyNb-z1TIyz1≤-2z1TK1z1-λminΓ1z2Tz2-STK2S-λminΘ1zfTzf+M-1¯z2zτ+-Γ1-I¯z2zf+N¯czszsS+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1-NbTIyNbz1TK1z1NbTIyNb-z1TIyz1To proceed to the following analysis, Young’s inequality is considered:(31)ab≤ι4a2+1ιb2where a and b denote vectors and ι represents positive constant. With this inequality, the N¯c(zs)zsS term in (30) becomes(32)N¯czszsS≤σ34N¯c2zszs2+1σ3S2≤σ34N¯c2zsz12+z22+S2+zτ2+zf2+1σ3S2Similar property is held for other terms in (30). In addition, considering the requirement of z1<Nb, the following inequalities yield:(33)σ38Nc¯2zsz1Tz1-z1TK1z1-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1≤-λminK1-σ3/8Nc¯2zsNbTIxNbz1TIxz1NbTIxNb-z1TIxz1≤-λminK1-σ38Nc¯2zsNbTIxNblog⁡NbTIxNbNbTIxNb-z1TIxz1For y, the similar transformation can be derived. Define(34)Ξx=-λminK1-σ38Nc¯2zsNbTIxNblog⁡NbTIxNbNbTIxNb-z1TIxz1Ξy=-λminK1-σ38Nc¯2zsNbTIyNblog⁡NbTIyNbNbTIyNb-z1TIyz1Invoking (31), (32), (33), and (34), (30) can be deducted as(35)V˙2≤Ξx+Ξy-σ34Nc¯2zsz2Tz2-λminK2-σ34Nc¯2zsSTS-λminΘ1-σ34Nc¯2zszfTzf+σ1M-1¯24z22+1σ1+σ34Nc¯2zszτ2+σ2-Γ1-I¯24z22+1σ2zf2+1σ3S2+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθCauchy-Schwarz inequality can help to regulate the upper bound of zτ as(36)zτ2≤τ1∫t-τ1tV˙g1θ2dθMoreover, the following fact can be proven:(37)∫t-τ1t∫wtV˙g1θ2dθdw≤τ1∫t-τ1tV˙g1θ2dθConsidering (36) and (37), (35) gives(38)V˙2≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf+υτ1V˙g1θ2-υ-τ1σ1-τ1σ34Nc¯2zs∫t-τ1tV˙g1θ2dθ≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf-υτ1-1σ1-σ34Nc¯2zs∫t-τ1t∫wtV˙g1θ2dθdw+υτ1V˙g1θ2≤-ρcV2+βcwhere ρc,βc>0 and they satisfy ρc=min[2Ξx,2Ξy,2(λmin(Γ1)-σ1M-1¯2/4-(σ3/4)Nc¯2(zs)-σ2(-Γ1-I)¯2/4)2(λmin(K2)-1/σ3 − (σ3/4)Nc¯2(zs))/λmax(M), 2(λmin(Θ1)-1/σ2-(σ3/4)Nc¯2(zs)), and (1/τ1-1/σ1υ-(σ3/4υ)Nc¯2(zs))] and the tuning parameters are selected λmin(K1)>(σ3/12)Nc¯2(zs), λmin(Γ1)+σ2-Γ1-I¯2/4>σ1M-1¯2/4+(σ3/4)Nc¯2(zs), λmin(K2)>1/σ3+(σ3/4)Nc¯2(zs), λmin(Θ1)>1/σ2+(σ3/4)Nc¯2(zs), and υ/τ1>1/σ1+(σ3/4)Nc¯2(zs). βc=υτ1V˙g1(θ)2.Lemma 8 (see [26, 27]).
Given bounded initial conditions, if there exists aC1 continuous and positive definite Lyapunov function V(x) satisfying v1(x)≤V(x)≤v2(x), such that V˙≤-αV(x)+β, where v1, v2:, Rn→R are class K functions and α,β>0, then the solution x(t) is uniformly bounded.Remark 9.
Combining Lemma8, Remark 7 and (13), (14), (15), (16), (17), (18), (19), (20), (21), (22), (23), (24), (25), (26), (28), (29), (30), (31), (32), (33), (34), (35), (36), (37), and (38), the Semiglobally Uniform Boundness (SGUB) of all the close-loop signals is validated with the involvement of input delay. Moreover, the tracking error can be constrained within z1≤Nb.


#### 3.2.3. Control Design with System Disturbance

Due to the system parameter uncertainty and external disturbance, the dynamic model is augmented with a disturbance term in the following format.(39)MTdg2,Tg2,tT˙dg2+CTdg2,Tg2,tTdg2+gTg2,t=Vg1t-τ1+dTdg2,Tg2,tUnder the disturbance case, the first step of the control design is totally similar to Step 1. And the whole verification before (23) remains the same; (24) will be adjusted as(40)MS˙=Mα˙c+Ms-d+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2Assumption 10.
The disturbanced(Tdg2,Tg2,t) considered in this paper is bounded such that a global approximator can be employed to estimate this signal.In order to approximate the unknown disturbance, a RBF neural network is employed.(41)d=W^cTScZc+ϵcDefine W^c,Wc∗,ϵc as the estimated weights, optimal weights, and approximation error, respectively. Zc stands for the input vector of the neural network. The specific expression on Zc will be introduced in the simulation section. Develop the update law of the NN weights as(42)W^˙ci=-ΥiSciZcSi+ξiW^ciControl input under this condition should be augmented into(43)Vg1bt=Vg1t-W^cTScZcThe control law in (43) will be guaranteed with SGUB of all the close-loop system signals.Proof.
Denote the optimal weights asW~c=W^c-Wc∗. The Lyapunov candidate function in the presence of disturbance can be selected as(44)V3=V2+12∑i=13W~ciTΥ-1W~ciConsidering the weight updating law (42), time derivative of V3 becomes(45)V˙3=V˙2-STd+STW^cScZc+∑i=12W~ciTΥ-1W~˙ci=V˙2-STd+STW^cScZc+∑i=12W~ciTΥ-1W^˙ci=V˙2-STd+STW^cScZc-∑i=12W~ciTSciZcSi-∑i=12ξiW~ciTW^ci=V˙2-STd+STW^cScZc-∑i=12W^ciTSciZcSi+∑i=12Wci∗TSciZcSi-∑i=12ξiW~ciTW^ci=V˙2-∑i=12ξiW~ciTW^ciIncorporating the inequality(46)-ξiW~ciTW^ci=-ξiW~ciTW~ci+Wci∗=-ξiW~ci2-ξiW^ci-Wci∗Wci∗=-ξiW~ci2-ξiW~ci+Wci∗Wci∗+ξiWci∗2=-ξi2W~ci2+ξi2Wci∗2-ξi2W~ci2+2W~ciWci∗+Wci∗2=-ξi2W~ci2+ξi2Wci∗2-ξi2W~ci+Wci∗2≤-ξi2W~ci2+ξi2Wci∗2V˙3 becomes(47)V˙3≤V˙2-∑i=12ξi2W~ci2+∑i=12ξi2Wci∗2≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf-υτ1-1σ1-σ34Nc¯2zs∫t-τ1t∫wtV˙g1θ2dθdw-12∑i=13ξiλmaxΥ-1W~ciTΥ-1W~ci+υτ1V˙g1θ2+∑i=13ξi2Wci∗2≤ρc′V3+βc′where ρc′,βc′>0 and they satisfy ρc′=min[2Ξx,2Ξy,2(λmin(Γ1)-σ1M-1¯2/4-(σ3/4)Nc¯2(zs)-σ2-Γ1-I¯2/4), 2(λmin(K2)-1/σ3-(σ3/4)Nc¯2(zs))/λmax(M), 2(λmin(Θ1)-1/σ2-(σ3/4)Nc¯2(zs)), and (1/τ1-1/σ1υ-(σ3/4υ)Nc¯2(zs))], the tuning parameters are selected λmin(K1)>(σ3/12)Nc¯2(zs), λmin(Γ1)+σ2-Γ1-I¯2/4>σ1M-1¯2/4+(σ3/4)Nc¯2(zs), λmin(K2)>1/σ3+(σ3/4)Nc¯2(zs), λmin(Θ1)>1/σ2+(σ3/4)Nc¯2(zs), and υ/τ1>1/σ1+(σ3/4)Nc¯2(zs). βc′=υτ1V˙g1(θ)2+σ5|ϵc¯∥2+∑i=12‍ξi/2Wci∗2.
With similar deduction of Remark9, the SGUB can still be guaranteed after the augmentation of NN control.


## 3.1. Preliminaries on NN Approximation

Radial Basis Function Neural Network (RBFNN) has been demonstrated to have outstanding function approximation and learning capabilities. If we have a continuous nonlinear functionϕ(z):Rn→Rm defined on a compact set Ωz, it can be approximated by RBFNN with following expression:(5)ϕNNZ=WTSNNZwhere Z∈Ωz⊂Rn represents the input vector. W∈Rp×m stands for the adjustable weight matrix and p is the number of neuron. SNN(Z)=[s1(Z),…,sp(Z)]T is the basis function vector, where si(Z)=exp⁡[-(Z-μi)T(Z-μi)/ψi2],i=1,2,…,p. μi=[μi1,μi2,…,μin] are the centers of respective field and ψi denotes Gaussian function width. Concretely, (5) can be described as(6)ϕZ=W∗TSNNZ+εZwhere W∗ is the ideal weight and ε(Z)∈Rm denotes the approximation error and is subject to ε(Z)≤σ, where σ should be a positive constant. The ideal weight matrix W∗ can be defined as(7)W∗≔arg⁡minW∈Rp×m⁡supZ∈ΩZ⁡ϕZ-WSNNZIt should be noted that W∗ is only an “artificial" quantity for the purpose of analysis. In control design, the ideal values will be estimated by W with developed updating law [20]. Moreover, NN approximation error ε(Z) indicates the minimum deviation between optimal approximation solution, i.e., W∗SNN(Z) and the unknown function ϕ(Z). It has been proved that approximation with neural networks can theoretically have any expected accuracy if the neuron number is large enough; i.e., ε(Z) is able to be arbitrarily small if p is sufficiently large [20].


## 3.2. Adaptive Control with Multiple Constraints

### 3.2.1. Mathematical Model Formulation

If we are already given the supply water temperatureTg2 in secondary network, our main aim is to control the valve opening degree in primary network to guarantee the supply water temperature following the desired values. Concretely, the given temperature value Tg2 acts as the desired signal symbol as Tg2d. The transfer function in this heat exchange process can be expressed as [22](8)G1s=KpT1s+1T2s+1e-τ1swhere Kp denotes a scalar coefficient, T1 and T2 represent the constant inertial time, and τ1 is the lagging time. This mathematical model depicts the transformation from system input Vg1, i.e., opening degree of the valve in primary pipe network, to the system output Tg2, i.e., actual water temperature in secondary pipe network. For the purpose of control design, we convert the model above into a time-domain formulation.(9)T1T2KpT¨g2+T1+T2KpT˙g2+Tg2Kp=Vg1t-τ1Equation (9) can be concisely written as(10)MT¨g2+CT˙g2+gTg2=Vg1t-τ1where M=T1T2/Kp, C=T1+T2/Kp, and g(Tg2)=Tg2/Kp. In order to extend the generalization of our control, we consider the following remarks and assumptions.Remark 3.
Even though the parametersKp, T1, and T2 in [22] are assumed as constants, in practical application, it is apparent that the heat exchange process is time-varying. In other words, the aforementioned coefficients should be state-dependent and time-varying. In this regard, M, C, and g(Tg2) can be more accurately considered as M(T˙g2,Tg2,t), C(T˙g2,Tg2,t), and g(Tg2,t).Remark 4.
In actual system, external disturbance as well as the model uncertainty should be considered in model (10). For simplicity purpose, they are merged into single term as d(Tdg2,Tg2,t).Besides the remarks and assumption made above, one more variableTdg2=T˙g2 is considered which is actually denotes the temperature variation in secondary loop. Thus, the system model of (10) can be modified as(11)T˙g2=Tdg2MTdg2,Tg2,tT˙dg2+CTdg2,Tg2,tTdg2+gTg2,t=Vg1t-τ1+dTdg2,Tg2,tAssumption 5.
In many researches, only single primary/secondary pipe network is considered; i.e., aforementioned model is defined as an SISO system. In this paper, system is extended to be an MIMO system which can produce multiple water temperatures in secondary pipe network to satisfy different operating requirements.Assumption 6.
The inertia matrixM is invertible and M-1 is bounded. The upper bound will be denoted as M-1∞≤M-1¯, where M-1¯∈R is a positive constant bound.With Assumption5, the input and output variables become a second order vector Tg2=[Tg21,Tg22]T and Vg1=[Vg11,Vg12]T. Additionally, since the disturbance term d(Tdg2,Tg2,t) is going to be estimated with an NN observer, it is neglected temporarily. The resultant system model can be concisely expressed as follows. Herein, we slightly abuse the use of M(Tdg2,Tg2,t), C(Tdg2,Tg2,t) to M, C(12)T˙g2=Tdg2MT˙dg2+CTdg2+g=Vg1t-τ1


### 3.2.2. Adaptive Control Design without Disturbance

In this paper, the output constraints, i.e., the tracked secondary loop water temperature errors, will also be considered to regulate within a desired range. In order to handle the output constraint problem, Symmetry Barrier Lyapunov Function (SBLF) [23] is employed together with a backstepping approach for the control development.Step 1.
Denote(13)z1=Tg2d-Tg2and(14)z2=αc-Tdg2where the desired trajectory satisfies Tg2d,T˙g2d∈L∞. αc refers to stabilizing function. Select a positive definite and C1 continuous SBLF candidate as(15)V1=12logNbTIxNbNbTIxNb-z1TIxz1+12logNbTIyNbNbTIyNb-z1TIyz1where(16)Ix=1000Iy=0001Nb∈R2×1 denotes the tracking error constraint such that z1≤Nb should be satisfied through the control.Remark 7.
In practical application, the initial conditions of water temperature and temperature variation are the same with the desired values. Hence,z10<Nb should be satisfied.Time derivative ofV1 gives(17)V˙1=z1TIxz˙1NbTIxNb-z1TIxz1+z1TIyz˙1NbTIyNb-z1TIyz1Differentiating z1 with respect to time gives(18)z˙1=T˙g2d-αc-z2Substituting (18) into (17), we obtain(19)V˙1=z1TIxT˙g2d-αc-z2NbTIxNb-z1TIxz1+z1TIyT˙g2d-αc-z2NbTIyNb-z1TIyz1Design the stabling function αc as(20)αc=T˙g2d+NbTNb-z1Tz1K1z1Substituteing (20) into (19), the following can be derived:(21)V˙1=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1Step 2.
Invoke an auxiliary stateS∈R2×1 to cope with the input delay which has the expression following formulation:(22)S=z2-M-1∫t-τ1tVg1θdθ-zfwhere zf∈R2×1 is developed with the following adaptive law.(23)z˙f=K2S-Γ1z2-Θzfwhere K2,Γ1,Θ∈R2×2 are positive coefficients. Multiplying both sides of (22) with M and involving Ms=C(Tdg2,Tg2,t)Tdg2+g(Tg2,t), the derivative of MS gives(24)MS˙=Mz˙2-Vg1t+Vg1t-τ1-z˙f=Mα˙c+CTdg2,Tg2,tTdg2+gTg2,t-Vg1t-K2S+Θzf+Γ1z2=Mα˙c+Ms+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2where without loss of generality M is considered as a constant matrix since the inertial time constants change very slowly. The time-varying part can be added to the disturbance term d. Moreover, Nc is designed as follows with the consideration of Mean Value Theorem [24]:(25)Nc=Θzf+Γ1z2+K2z2+ST+S˙Tz2Nc≤N¯czszswhere the bounding function N¯c(zs) denotes a globally positive function. zs is defined as zs=[z1T,z2T,ST,zτT,zfT]T, where zτ∈R2×1 refers to(26)zτ=Vg1t-Vg1t-τ1=∫t-τ1tV˙g1θdθAfter the introduction of auxiliary state S, the input-delay system is transferred into a delay-free system as expressed in (24). In order to prove the stability of the resultant system, a quadratic form Lyapunov-Krasovskii candidate function can be designated as [25](27)V2=V1+12z2Tz2+12STMS+12zfTzf+υ∫t-τ1t∫wtV˙g1θ2dθdw
DifferentiateV2 and incorporate (21), (22), (23), and (24); this yields(28)V˙2=V˙1+z2Tz˙2+STMS˙+zfz˙f+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1+z2TS˙-M-1Vg1t-τ1-Vg1t+K2S-Θzf-Γ1z2+STMα˙c+Ms+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2+zfTK2S-zfTΘzf-zfTΓ1z2+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1-z2TΓ1z2+z2TM-1zτ-z2TΓ1+Izf-STK2S+STMα˙c+Ms+Nc-Vg1θ+zfTK2S-zfTΘzf+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθThe control law can be developed as follows:(29)Vg1t=Mα˙c+Ms+K2zf+ST+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1+NbTIxNbz1TK1z1NbTIxNb-z1TIxz1+NbTIyNbz1TK1z1NbTIyNb-z1TIyz1Substituting (29) into (28) and invoking (25) and Assumption 6, we have(30)V˙2=-2z1TK1z1-z2TΓ1z2-STK2S-zfTΘzf+z2TM-1zτ-z2TΓ1+Izf+STNc+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1-NbTIyNbz1TK1z1NbTIyNb-z1TIyz1≤-2z1TK1z1-λminΓ1z2Tz2-STK2S-λminΘ1zfTzf+M-1¯z2zτ+-Γ1-I¯z2zf+N¯czszsS+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1-NbTIyNbz1TK1z1NbTIyNb-z1TIyz1To proceed to the following analysis, Young’s inequality is considered:(31)ab≤ι4a2+1ιb2where a and b denote vectors and ι represents positive constant. With this inequality, the N¯c(zs)zsS term in (30) becomes(32)N¯czszsS≤σ34N¯c2zszs2+1σ3S2≤σ34N¯c2zsz12+z22+S2+zτ2+zf2+1σ3S2Similar property is held for other terms in (30). In addition, considering the requirement of z1<Nb, the following inequalities yield:(33)σ38Nc¯2zsz1Tz1-z1TK1z1-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1≤-λminK1-σ3/8Nc¯2zsNbTIxNbz1TIxz1NbTIxNb-z1TIxz1≤-λminK1-σ38Nc¯2zsNbTIxNblog⁡NbTIxNbNbTIxNb-z1TIxz1For y, the similar transformation can be derived. Define(34)Ξx=-λminK1-σ38Nc¯2zsNbTIxNblog⁡NbTIxNbNbTIxNb-z1TIxz1Ξy=-λminK1-σ38Nc¯2zsNbTIyNblog⁡NbTIyNbNbTIyNb-z1TIyz1Invoking (31), (32), (33), and (34), (30) can be deducted as(35)V˙2≤Ξx+Ξy-σ34Nc¯2zsz2Tz2-λminK2-σ34Nc¯2zsSTS-λminΘ1-σ34Nc¯2zszfTzf+σ1M-1¯24z22+1σ1+σ34Nc¯2zszτ2+σ2-Γ1-I¯24z22+1σ2zf2+1σ3S2+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθCauchy-Schwarz inequality can help to regulate the upper bound of zτ as(36)zτ2≤τ1∫t-τ1tV˙g1θ2dθMoreover, the following fact can be proven:(37)∫t-τ1t∫wtV˙g1θ2dθdw≤τ1∫t-τ1tV˙g1θ2dθConsidering (36) and (37), (35) gives(38)V˙2≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf+υτ1V˙g1θ2-υ-τ1σ1-τ1σ34Nc¯2zs∫t-τ1tV˙g1θ2dθ≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf-υτ1-1σ1-σ34Nc¯2zs∫t-τ1t∫wtV˙g1θ2dθdw+υτ1V˙g1θ2≤-ρcV2+βcwhere ρc,βc>0 and they satisfy ρc=min[2Ξx,2Ξy,2(λmin(Γ1)-σ1M-1¯2/4-(σ3/4)Nc¯2(zs)-σ2(-Γ1-I)¯2/4)2(λmin(K2)-1/σ3 − (σ3/4)Nc¯2(zs))/λmax(M), 2(λmin(Θ1)-1/σ2-(σ3/4)Nc¯2(zs)), and (1/τ1-1/σ1υ-(σ3/4υ)Nc¯2(zs))] and the tuning parameters are selected λmin(K1)>(σ3/12)Nc¯2(zs), λmin(Γ1)+σ2-Γ1-I¯2/4>σ1M-1¯2/4+(σ3/4)Nc¯2(zs), λmin(K2)>1/σ3+(σ3/4)Nc¯2(zs), λmin(Θ1)>1/σ2+(σ3/4)Nc¯2(zs), and υ/τ1>1/σ1+(σ3/4)Nc¯2(zs). βc=υτ1V˙g1(θ)2.Lemma 8 (see [26, 27]).
Given bounded initial conditions, if there exists aC1 continuous and positive definite Lyapunov function V(x) satisfying v1(x)≤V(x)≤v2(x), such that V˙≤-αV(x)+β, where v1, v2:, Rn→R are class K functions and α,β>0, then the solution x(t) is uniformly bounded.Remark 9.
Combining Lemma8, Remark 7 and (13), (14), (15), (16), (17), (18), (19), (20), (21), (22), (23), (24), (25), (26), (28), (29), (30), (31), (32), (33), (34), (35), (36), (37), and (38), the Semiglobally Uniform Boundness (SGUB) of all the close-loop signals is validated with the involvement of input delay. Moreover, the tracking error can be constrained within z1≤Nb.


### 3.2.3. Control Design with System Disturbance

Due to the system parameter uncertainty and external disturbance, the dynamic model is augmented with a disturbance term in the following format.(39)MTdg2,Tg2,tT˙dg2+CTdg2,Tg2,tTdg2+gTg2,t=Vg1t-τ1+dTdg2,Tg2,tUnder the disturbance case, the first step of the control design is totally similar to Step 1. And the whole verification before (23) remains the same; (24) will be adjusted as(40)MS˙=Mα˙c+Ms-d+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2Assumption 10.
The disturbanced(Tdg2,Tg2,t) considered in this paper is bounded such that a global approximator can be employed to estimate this signal.In order to approximate the unknown disturbance, a RBF neural network is employed.(41)d=W^cTScZc+ϵcDefine W^c,Wc∗,ϵc as the estimated weights, optimal weights, and approximation error, respectively. Zc stands for the input vector of the neural network. The specific expression on Zc will be introduced in the simulation section. Develop the update law of the NN weights as(42)W^˙ci=-ΥiSciZcSi+ξiW^ciControl input under this condition should be augmented into(43)Vg1bt=Vg1t-W^cTScZcThe control law in (43) will be guaranteed with SGUB of all the close-loop system signals.Proof.
Denote the optimal weights asW~c=W^c-Wc∗. The Lyapunov candidate function in the presence of disturbance can be selected as(44)V3=V2+12∑i=13W~ciTΥ-1W~ciConsidering the weight updating law (42), time derivative of V3 becomes(45)V˙3=V˙2-STd+STW^cScZc+∑i=12W~ciTΥ-1W~˙ci=V˙2-STd+STW^cScZc+∑i=12W~ciTΥ-1W^˙ci=V˙2-STd+STW^cScZc-∑i=12W~ciTSciZcSi-∑i=12ξiW~ciTW^ci=V˙2-STd+STW^cScZc-∑i=12W^ciTSciZcSi+∑i=12Wci∗TSciZcSi-∑i=12ξiW~ciTW^ci=V˙2-∑i=12ξiW~ciTW^ciIncorporating the inequality(46)-ξiW~ciTW^ci=-ξiW~ciTW~ci+Wci∗=-ξiW~ci2-ξiW^ci-Wci∗Wci∗=-ξiW~ci2-ξiW~ci+Wci∗Wci∗+ξiWci∗2=-ξi2W~ci2+ξi2Wci∗2-ξi2W~ci2+2W~ciWci∗+Wci∗2=-ξi2W~ci2+ξi2Wci∗2-ξi2W~ci+Wci∗2≤-ξi2W~ci2+ξi2Wci∗2V˙3 becomes(47)V˙3≤V˙2-∑i=12ξi2W~ci2+∑i=12ξi2Wci∗2≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf-υτ1-1σ1-σ34Nc¯2zs∫t-τ1t∫wtV˙g1θ2dθdw-12∑i=13ξiλmaxΥ-1W~ciTΥ-1W~ci+υτ1V˙g1θ2+∑i=13ξi2Wci∗2≤ρc′V3+βc′where ρc′,βc′>0 and they satisfy ρc′=min[2Ξx,2Ξy,2(λmin(Γ1)-σ1M-1¯2/4-(σ3/4)Nc¯2(zs)-σ2-Γ1-I¯2/4), 2(λmin(K2)-1/σ3-(σ3/4)Nc¯2(zs))/λmax(M), 2(λmin(Θ1)-1/σ2-(σ3/4)Nc¯2(zs)), and (1/τ1-1/σ1υ-(σ3/4υ)Nc¯2(zs))], the tuning parameters are selected λmin(K1)>(σ3/12)Nc¯2(zs), λmin(Γ1)+σ2-Γ1-I¯2/4>σ1M-1¯2/4+(σ3/4)Nc¯2(zs), λmin(K2)>1/σ3+(σ3/4)Nc¯2(zs), λmin(Θ1)>1/σ2+(σ3/4)Nc¯2(zs), and υ/τ1>1/σ1+(σ3/4)Nc¯2(zs). βc′=υτ1V˙g1(θ)2+σ5|ϵc¯∥2+∑i=12‍ξi/2Wci∗2.
With similar deduction of Remark9, the SGUB can still be guaranteed after the augmentation of NN control.


## 3.2.1. Mathematical Model Formulation

If we are already given the supply water temperatureTg2 in secondary network, our main aim is to control the valve opening degree in primary network to guarantee the supply water temperature following the desired values. Concretely, the given temperature value Tg2 acts as the desired signal symbol as Tg2d. The transfer function in this heat exchange process can be expressed as [22](8)G1s=KpT1s+1T2s+1e-τ1swhere Kp denotes a scalar coefficient, T1 and T2 represent the constant inertial time, and τ1 is the lagging time. This mathematical model depicts the transformation from system input Vg1, i.e., opening degree of the valve in primary pipe network, to the system output Tg2, i.e., actual water temperature in secondary pipe network. For the purpose of control design, we convert the model above into a time-domain formulation.(9)T1T2KpT¨g2+T1+T2KpT˙g2+Tg2Kp=Vg1t-τ1Equation (9) can be concisely written as(10)MT¨g2+CT˙g2+gTg2=Vg1t-τ1where M=T1T2/Kp, C=T1+T2/Kp, and g(Tg2)=Tg2/Kp. In order to extend the generalization of our control, we consider the following remarks and assumptions.Remark 3.
Even though the parametersKp, T1, and T2 in [22] are assumed as constants, in practical application, it is apparent that the heat exchange process is time-varying. In other words, the aforementioned coefficients should be state-dependent and time-varying. In this regard, M, C, and g(Tg2) can be more accurately considered as M(T˙g2,Tg2,t), C(T˙g2,Tg2,t), and g(Tg2,t).Remark 4.
In actual system, external disturbance as well as the model uncertainty should be considered in model (10). For simplicity purpose, they are merged into single term as d(Tdg2,Tg2,t).Besides the remarks and assumption made above, one more variableTdg2=T˙g2 is considered which is actually denotes the temperature variation in secondary loop. Thus, the system model of (10) can be modified as(11)T˙g2=Tdg2MTdg2,Tg2,tT˙dg2+CTdg2,Tg2,tTdg2+gTg2,t=Vg1t-τ1+dTdg2,Tg2,tAssumption 5.
In many researches, only single primary/secondary pipe network is considered; i.e., aforementioned model is defined as an SISO system. In this paper, system is extended to be an MIMO system which can produce multiple water temperatures in secondary pipe network to satisfy different operating requirements.Assumption 6.
The inertia matrixM is invertible and M-1 is bounded. The upper bound will be denoted as M-1∞≤M-1¯, where M-1¯∈R is a positive constant bound.With Assumption5, the input and output variables become a second order vector Tg2=[Tg21,Tg22]T and Vg1=[Vg11,Vg12]T. Additionally, since the disturbance term d(Tdg2,Tg2,t) is going to be estimated with an NN observer, it is neglected temporarily. The resultant system model can be concisely expressed as follows. Herein, we slightly abuse the use of M(Tdg2,Tg2,t), C(Tdg2,Tg2,t) to M, C(12)T˙g2=Tdg2MT˙dg2+CTdg2+g=Vg1t-τ1


## 3.2.2. Adaptive Control Design without Disturbance

In this paper, the output constraints, i.e., the tracked secondary loop water temperature errors, will also be considered to regulate within a desired range. In order to handle the output constraint problem, Symmetry Barrier Lyapunov Function (SBLF) [23] is employed together with a backstepping approach for the control development.Step 1.
Denote(13)z1=Tg2d-Tg2and(14)z2=αc-Tdg2where the desired trajectory satisfies Tg2d,T˙g2d∈L∞. αc refers to stabilizing function. Select a positive definite and C1 continuous SBLF candidate as(15)V1=12logNbTIxNbNbTIxNb-z1TIxz1+12logNbTIyNbNbTIyNb-z1TIyz1where(16)Ix=1000Iy=0001Nb∈R2×1 denotes the tracking error constraint such that z1≤Nb should be satisfied through the control.Remark 7.
In practical application, the initial conditions of water temperature and temperature variation are the same with the desired values. Hence,z10<Nb should be satisfied.Time derivative ofV1 gives(17)V˙1=z1TIxz˙1NbTIxNb-z1TIxz1+z1TIyz˙1NbTIyNb-z1TIyz1Differentiating z1 with respect to time gives(18)z˙1=T˙g2d-αc-z2Substituting (18) into (17), we obtain(19)V˙1=z1TIxT˙g2d-αc-z2NbTIxNb-z1TIxz1+z1TIyT˙g2d-αc-z2NbTIyNb-z1TIyz1Design the stabling function αc as(20)αc=T˙g2d+NbTNb-z1Tz1K1z1Substituteing (20) into (19), the following can be derived:(21)V˙1=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1Step 2.
Invoke an auxiliary stateS∈R2×1 to cope with the input delay which has the expression following formulation:(22)S=z2-M-1∫t-τ1tVg1θdθ-zfwhere zf∈R2×1 is developed with the following adaptive law.(23)z˙f=K2S-Γ1z2-Θzfwhere K2,Γ1,Θ∈R2×2 are positive coefficients. Multiplying both sides of (22) with M and involving Ms=C(Tdg2,Tg2,t)Tdg2+g(Tg2,t), the derivative of MS gives(24)MS˙=Mz˙2-Vg1t+Vg1t-τ1-z˙f=Mα˙c+CTdg2,Tg2,tTdg2+gTg2,t-Vg1t-K2S+Θzf+Γ1z2=Mα˙c+Ms+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2where without loss of generality M is considered as a constant matrix since the inertial time constants change very slowly. The time-varying part can be added to the disturbance term d. Moreover, Nc is designed as follows with the consideration of Mean Value Theorem [24]:(25)Nc=Θzf+Γ1z2+K2z2+ST+S˙Tz2Nc≤N¯czszswhere the bounding function N¯c(zs) denotes a globally positive function. zs is defined as zs=[z1T,z2T,ST,zτT,zfT]T, where zτ∈R2×1 refers to(26)zτ=Vg1t-Vg1t-τ1=∫t-τ1tV˙g1θdθAfter the introduction of auxiliary state S, the input-delay system is transferred into a delay-free system as expressed in (24). In order to prove the stability of the resultant system, a quadratic form Lyapunov-Krasovskii candidate function can be designated as [25](27)V2=V1+12z2Tz2+12STMS+12zfTzf+υ∫t-τ1t∫wtV˙g1θ2dθdw
DifferentiateV2 and incorporate (21), (22), (23), and (24); this yields(28)V˙2=V˙1+z2Tz˙2+STMS˙+zfz˙f+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1+z2TS˙-M-1Vg1t-τ1-Vg1t+K2S-Θzf-Γ1z2+STMα˙c+Ms+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2+zfTK2S-zfTΘzf-zfTΓ1z2+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ=-2z1TK1z1+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1-z2TΓ1z2+z2TM-1zτ-z2TΓ1+Izf-STK2S+STMα˙c+Ms+Nc-Vg1θ+zfTK2S-zfTΘzf+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθThe control law can be developed as follows:(29)Vg1t=Mα˙c+Ms+K2zf+ST+z1TIxz2NbTIxNb-z1TIxz1+z1TIyz2NbTIyNb-z1TIyz1+NbTIxNbz1TK1z1NbTIxNb-z1TIxz1+NbTIyNbz1TK1z1NbTIyNb-z1TIyz1Substituting (29) into (28) and invoking (25) and Assumption 6, we have(30)V˙2=-2z1TK1z1-z2TΓ1z2-STK2S-zfTΘzf+z2TM-1zτ-z2TΓ1+Izf+STNc+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1-NbTIyNbz1TK1z1NbTIyNb-z1TIyz1≤-2z1TK1z1-λminΓ1z2Tz2-STK2S-λminΘ1zfTzf+M-1¯z2zτ+-Γ1-I¯z2zf+N¯czszsS+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθ-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1-NbTIyNbz1TK1z1NbTIyNb-z1TIyz1To proceed to the following analysis, Young’s inequality is considered:(31)ab≤ι4a2+1ιb2where a and b denote vectors and ι represents positive constant. With this inequality, the N¯c(zs)zsS term in (30) becomes(32)N¯czszsS≤σ34N¯c2zszs2+1σ3S2≤σ34N¯c2zsz12+z22+S2+zτ2+zf2+1σ3S2Similar property is held for other terms in (30). In addition, considering the requirement of z1<Nb, the following inequalities yield:(33)σ38Nc¯2zsz1Tz1-z1TK1z1-NbTIxNbz1TK1z1NbTIxNb-z1TIxz1≤-λminK1-σ3/8Nc¯2zsNbTIxNbz1TIxz1NbTIxNb-z1TIxz1≤-λminK1-σ38Nc¯2zsNbTIxNblog⁡NbTIxNbNbTIxNb-z1TIxz1For y, the similar transformation can be derived. Define(34)Ξx=-λminK1-σ38Nc¯2zsNbTIxNblog⁡NbTIxNbNbTIxNb-z1TIxz1Ξy=-λminK1-σ38Nc¯2zsNbTIyNblog⁡NbTIyNbNbTIyNb-z1TIyz1Invoking (31), (32), (33), and (34), (30) can be deducted as(35)V˙2≤Ξx+Ξy-σ34Nc¯2zsz2Tz2-λminK2-σ34Nc¯2zsSTS-λminΘ1-σ34Nc¯2zszfTzf+σ1M-1¯24z22+1σ1+σ34Nc¯2zszτ2+σ2-Γ1-I¯24z22+1σ2zf2+1σ3S2+υτ1V˙g1θ2-υ∫t-τ1tV˙g1θ2dθCauchy-Schwarz inequality can help to regulate the upper bound of zτ as(36)zτ2≤τ1∫t-τ1tV˙g1θ2dθMoreover, the following fact can be proven:(37)∫t-τ1t∫wtV˙g1θ2dθdw≤τ1∫t-τ1tV˙g1θ2dθConsidering (36) and (37), (35) gives(38)V˙2≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf+υτ1V˙g1θ2-υ-τ1σ1-τ1σ34Nc¯2zs∫t-τ1tV˙g1θ2dθ≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf-υτ1-1σ1-σ34Nc¯2zs∫t-τ1t∫wtV˙g1θ2dθdw+υτ1V˙g1θ2≤-ρcV2+βcwhere ρc,βc>0 and they satisfy ρc=min[2Ξx,2Ξy,2(λmin(Γ1)-σ1M-1¯2/4-(σ3/4)Nc¯2(zs)-σ2(-Γ1-I)¯2/4)2(λmin(K2)-1/σ3 − (σ3/4)Nc¯2(zs))/λmax(M), 2(λmin(Θ1)-1/σ2-(σ3/4)Nc¯2(zs)), and (1/τ1-1/σ1υ-(σ3/4υ)Nc¯2(zs))] and the tuning parameters are selected λmin(K1)>(σ3/12)Nc¯2(zs), λmin(Γ1)+σ2-Γ1-I¯2/4>σ1M-1¯2/4+(σ3/4)Nc¯2(zs), λmin(K2)>1/σ3+(σ3/4)Nc¯2(zs), λmin(Θ1)>1/σ2+(σ3/4)Nc¯2(zs), and υ/τ1>1/σ1+(σ3/4)Nc¯2(zs). βc=υτ1V˙g1(θ)2.Lemma 8 (see [26, 27]).
Given bounded initial conditions, if there exists aC1 continuous and positive definite Lyapunov function V(x) satisfying v1(x)≤V(x)≤v2(x), such that V˙≤-αV(x)+β, where v1, v2:, Rn→R are class K functions and α,β>0, then the solution x(t) is uniformly bounded.Remark 9.
Combining Lemma8, Remark 7 and (13), (14), (15), (16), (17), (18), (19), (20), (21), (22), (23), (24), (25), (26), (28), (29), (30), (31), (32), (33), (34), (35), (36), (37), and (38), the Semiglobally Uniform Boundness (SGUB) of all the close-loop signals is validated with the involvement of input delay. Moreover, the tracking error can be constrained within z1≤Nb.


## 3.2.3. Control Design with System Disturbance

Due to the system parameter uncertainty and external disturbance, the dynamic model is augmented with a disturbance term in the following format.(39)MTdg2,Tg2,tT˙dg2+CTdg2,Tg2,tTdg2+gTg2,t=Vg1t-τ1+dTdg2,Tg2,tUnder the disturbance case, the first step of the control design is totally similar to Step 1. And the whole verification before (23) remains the same; (24) will be adjusted as(40)MS˙=Mα˙c+Ms-d+Nc-Vg1t-K2S-K2z2-ST+S˙Tz2Assumption 10.
The disturbanced(Tdg2,Tg2,t) considered in this paper is bounded such that a global approximator can be employed to estimate this signal.In order to approximate the unknown disturbance, a RBF neural network is employed.(41)d=W^cTScZc+ϵcDefine W^c,Wc∗,ϵc as the estimated weights, optimal weights, and approximation error, respectively. Zc stands for the input vector of the neural network. The specific expression on Zc will be introduced in the simulation section. Develop the update law of the NN weights as(42)W^˙ci=-ΥiSciZcSi+ξiW^ciControl input under this condition should be augmented into(43)Vg1bt=Vg1t-W^cTScZcThe control law in (43) will be guaranteed with SGUB of all the close-loop system signals.Proof.
Denote the optimal weights asW~c=W^c-Wc∗. The Lyapunov candidate function in the presence of disturbance can be selected as(44)V3=V2+12∑i=13W~ciTΥ-1W~ciConsidering the weight updating law (42), time derivative of V3 becomes(45)V˙3=V˙2-STd+STW^cScZc+∑i=12W~ciTΥ-1W~˙ci=V˙2-STd+STW^cScZc+∑i=12W~ciTΥ-1W^˙ci=V˙2-STd+STW^cScZc-∑i=12W~ciTSciZcSi-∑i=12ξiW~ciTW^ci=V˙2-STd+STW^cScZc-∑i=12W^ciTSciZcSi+∑i=12Wci∗TSciZcSi-∑i=12ξiW~ciTW^ci=V˙2-∑i=12ξiW~ciTW^ciIncorporating the inequality(46)-ξiW~ciTW^ci=-ξiW~ciTW~ci+Wci∗=-ξiW~ci2-ξiW^ci-Wci∗Wci∗=-ξiW~ci2-ξiW~ci+Wci∗Wci∗+ξiWci∗2=-ξi2W~ci2+ξi2Wci∗2-ξi2W~ci2+2W~ciWci∗+Wci∗2=-ξi2W~ci2+ξi2Wci∗2-ξi2W~ci+Wci∗2≤-ξi2W~ci2+ξi2Wci∗2V˙3 becomes(47)V˙3≤V˙2-∑i=12ξi2W~ci2+∑i=12ξi2Wci∗2≤Ξx+Ξy-λminΓ1-σ1M-1¯24-σ34Nc¯2zs-σ2-Γ1-I¯24z2Tz2-λminK2-1σ3-σ34Nc¯2zsSTS-λminΘ1-1σ2-σ34Nc¯2zszfTzf-υτ1-1σ1-σ34Nc¯2zs∫t-τ1t∫wtV˙g1θ2dθdw-12∑i=13ξiλmaxΥ-1W~ciTΥ-1W~ci+υτ1V˙g1θ2+∑i=13ξi2Wci∗2≤ρc′V3+βc′where ρc′,βc′>0 and they satisfy ρc′=min[2Ξx,2Ξy,2(λmin(Γ1)-σ1M-1¯2/4-(σ3/4)Nc¯2(zs)-σ2-Γ1-I¯2/4), 2(λmin(K2)-1/σ3-(σ3/4)Nc¯2(zs))/λmax(M), 2(λmin(Θ1)-1/σ2-(σ3/4)Nc¯2(zs)), and (1/τ1-1/σ1υ-(σ3/4υ)Nc¯2(zs))], the tuning parameters are selected λmin(K1)>(σ3/12)Nc¯2(zs), λmin(Γ1)+σ2-Γ1-I¯2/4>σ1M-1¯2/4+(σ3/4)Nc¯2(zs), λmin(K2)>1/σ3+(σ3/4)Nc¯2(zs), λmin(Θ1)>1/σ2+(σ3/4)Nc¯2(zs), and υ/τ1>1/σ1+(σ3/4)Nc¯2(zs). βc′=υτ1V˙g1(θ)2+σ5|ϵc¯∥2+∑i=12‍ξi/2Wci∗2.
With similar deduction of Remark9, the SGUB can still be guaranteed after the augmentation of NN control.


## 4. Simulation Study

In this simulation, a numerical case study is conducted on a dual subnetwork heating system in the presence of output constraints and input delay.


### 4.1. Constrained Control with Input Delay

In this subsection, the control performance of proposed control is investigated first. As mentioned above, the complicated nonlinear system model is employed and the mass matrixM is separated and the major constant part is considered explicitly in the model. The specific details of the system are considered as M=diag⁡([249.27,235.28]), C=diag⁡([26.24∗Tdg11,23.26∗Tdg22]), g=0.8197∗Tg2, and the disturbance and model uncertainty error d=[0.05∗Tdg21∗cos⁡0.5∗t,0.05∗Tdg22∗cos⁡(0.5∗t+π/2)]. The whole simulation time is set as 40s and the desired temperature trajectory of secondary network Tg2d is defined as [22.49+(20.3-22.49)∗sin⁡(π/2∗40t+π/2),22.8+(21.7-22.8)∗sin⁡(π/(2∗40)t+π/2)]T. This design tracking trajectory describes the increment of temperature in secondary network from 20.3/21.7°C to 22.49/22.8°C in 40s, respectively. The constant input delay τ1 is fixed at 2.1s which is quite challenging for the control system. The control output temperature constraints are Nb=[0.15,0.12]T°C. It should be noted that the output constraints can be selected in an arbitrary fashion; however, too narrow constraints will cause the system instability.As for the control parameters, the main coefficients ofK1, K2, Γ1, and Θ are selected as diag⁡([0.0594,0.0864]), diag⁡([0.011,0.008]), diag⁡([0.0232,0.0232]), diag⁡([0.001,0.001]), respectively. In the neural networks, the learning rate and coefficients are Υ1 = Υ2 = 0.1, ξ1 = ξ2 = 0.01. After trial and error, the number of neural nodes are set as 32. The input Zc includes Tg2, Tdg2, αc, and α˙c. The corresponding centres are distributed within [-5,5], [-2,2], [-1,1], and [-2,2]. The control performance is presented in Figures 3–5.Figure 3
Tracking control result ofTg2 using proposed control.
(a)
Temperature tracked ofTg21 in subnetwork 1 (b)
Temperature tracked ofTg22 in subnetwork 2Figure 4
Constrained control errors with output limitation.
(a)
Tracking error of subnetwork 1 (b)
Tracking error of subnetwork 2Figure 5
Control input and temperature variation speed.
(a)
Control input of proposed control (b)
Temptation variation speedFrom the result, it can be observed that the overall control performance is excellent with proposed control in the presence of output constraints and input delay as well as disturbance and model uncertainty. In Figure3, the control trajectory can follow the desired temperature trajectory well. In the tail few seconds of the simulation, the control error increases become during this period and the increment of desired temperature becomes faster. Even this, the proposed control can still successfully regulate the tracking error within the desired constraints as shown in Figure 4. The error goes up in Figure 4(a) and before crossing the constraint limit, the control regulates it to reduce the error still within the desired area, which verifies the efficiency of proposed control. Finally, Figure 5 shows the control input and the temperature variation speed, i.e., Tdg2. The result reveals that all the close-loop signals are bounded under the proposed control scheme.


### 4.2. Control Performance with PID

To further illustrate the advantages of proposed control, a comparative study using PID control is conducted under the same system settings with input delay and disturbance. Due to the involvement of these challenging factors, the stability of the PID control can hardly be achieved. After many tuning trials, the best performance that PID control can achieve is as Figures6 and 7 show. The corresponding PID parameters are Kp=diag⁡([9.6,0;0,8]), Kd=diag⁡([4,0;0,2]), and Ki=diag⁡([0.05,0;0,0.05]).Figure 6
PID control output.
(a)
PID tracking control of subnetwork1 (b)
PID tracking control of subnetwork1Figure 7
PID control tracking error and input.
(a)
PID tracking error (b)
PID control inputAlthough all the close-loop signals including control error and control input are bounded, the errors are very large compared to the proposed control. It demonstrates that, under such challenging situation, traditional simple control cannot provide satisfactory performance and the proposed control manages to handle the issues efficiently.


## 4.1. Constrained Control with Input Delay

In this subsection, the control performance of proposed control is investigated first. As mentioned above, the complicated nonlinear system model is employed and the mass matrixM is separated and the major constant part is considered explicitly in the model. The specific details of the system are considered as M=diag⁡([249.27,235.28]), C=diag⁡([26.24∗Tdg11,23.26∗Tdg22]), g=0.8197∗Tg2, and the disturbance and model uncertainty error d=[0.05∗Tdg21∗cos⁡0.5∗t,0.05∗Tdg22∗cos⁡(0.5∗t+π/2)]. The whole simulation time is set as 40s and the desired temperature trajectory of secondary network Tg2d is defined as [22.49+(20.3-22.49)∗sin⁡(π/2∗40t+π/2),22.8+(21.7-22.8)∗sin⁡(π/(2∗40)t+π/2)]T. This design tracking trajectory describes the increment of temperature in secondary network from 20.3/21.7°C to 22.49/22.8°C in 40s, respectively. The constant input delay τ1 is fixed at 2.1s which is quite challenging for the control system. The control output temperature constraints are Nb=[0.15,0.12]T°C. It should be noted that the output constraints can be selected in an arbitrary fashion; however, too narrow constraints will cause the system instability.As for the control parameters, the main coefficients ofK1, K2, Γ1, and Θ are selected as diag⁡([0.0594,0.0864]), diag⁡([0.011,0.008]), diag⁡([0.0232,0.0232]), diag⁡([0.001,0.001]), respectively. In the neural networks, the learning rate and coefficients are Υ1 = Υ2 = 0.1, ξ1 = ξ2 = 0.01. After trial and error, the number of neural nodes are set as 32. The input Zc includes Tg2, Tdg2, αc, and α˙c. The corresponding centres are distributed within [-5,5], [-2,2], [-1,1], and [-2,2]. The control performance is presented in Figures 3–5.Figure 3
Tracking control result ofTg2 using proposed control.
(a)
Temperature tracked ofTg21 in subnetwork 1 (b)
Temperature tracked ofTg22 in subnetwork 2Figure 4
Constrained control errors with output limitation.
(a)
Tracking error of subnetwork 1 (b)
Tracking error of subnetwork 2Figure 5
Control input and temperature variation speed.
(a)
Control input of proposed control (b)
Temptation variation speedFrom the result, it can be observed that the overall control performance is excellent with proposed control in the presence of output constraints and input delay as well as disturbance and model uncertainty. In Figure3, the control trajectory can follow the desired temperature trajectory well. In the tail few seconds of the simulation, the control error increases become during this period and the increment of desired temperature becomes faster. Even this, the proposed control can still successfully regulate the tracking error within the desired constraints as shown in Figure 4. The error goes up in Figure 4(a) and before crossing the constraint limit, the control regulates it to reduce the error still within the desired area, which verifies the efficiency of proposed control. Finally, Figure 5 shows the control input and the temperature variation speed, i.e., Tdg2. The result reveals that all the close-loop signals are bounded under the proposed control scheme.


## 4.2. Control Performance with PID

To further illustrate the advantages of proposed control, a comparative study using PID control is conducted under the same system settings with input delay and disturbance. Due to the involvement of these challenging factors, the stability of the PID control can hardly be achieved. After many tuning trials, the best performance that PID control can achieve is as Figures6 and 7 show. The corresponding PID parameters are Kp=diag⁡([9.6,0;0,8]), Kd=diag⁡([4,0;0,2]), and Ki=diag⁡([0.05,0;0,0.05]).Figure 6
PID control output.
(a)
PID tracking control of subnetwork1 (b)
PID tracking control of subnetwork1Figure 7
PID control tracking error and input.
(a)
PID tracking error (b)
PID control inputAlthough all the close-loop signals including control error and control input are bounded, the errors are very large compared to the proposed control. It demonstrates that, under such challenging situation, traditional simple control cannot provide satisfactory performance and the proposed control manages to handle the issues efficiently.


## 5. Conclusion

In this paper, a novel adaptive NN based constrained control has been proposed to a city central heating system in the presence of input delay. The valve opening degree in the primary loop is controlled such that the supply water temperature in secondary loop can follow the desired values. In order to handle output tracking error constraints problem, BLF method is involved with designed Lyapunov candidate function. The input delay due to the slow action of valves and system signal transmission is resolved with an auxiliary system which facilitates to convert the delay system into a delay-free one. Finally, the unknown external disturbance and model uncertainty are approximated through neural network with developed adaptation law. Thorough simulation study is carried out to demonstrate that the proposed control can outperform traditional control scheme in such operation conditions.



---

*Source: 1714961-2018-11-19.xml*