Model instalacji CO w budynku B1

 

 

 

 

Oznaczenia:

 

Tzco1 – temperatura wody zasilającej halę technologiczną przy budynku B1, [°C],

Tzco2 – temperatura wody zasilającej budynek B1, [°C],

Tpco1 – temperatura wody powrotnej z hali technologicznej przy budynku B1, [°C],

Tpco2 – temperatura wody powrotnej z budynku B1, [°C],

Tpom1 – średnia temperatura w pomieszczeniach hali technologicznej, [°C],

Tpom2 – średnia temperatura w pomieszczeniach w budynku B1, [°C],

Fzco1 – przepływ wody zasilającej halę technologiczną, [t/h],

Fzco2 – przepływ wody zasilającej budynek B1, [t/h],

Uzco1 – stopień otwarcia zaworu regulującego przepływ wody zasilającej halę technologiczną, [0-100%],

Uzco2 – stopień otwarcia zaworu regulującego przepływ wody zasilającej budynek B1, [0-100%].

cw – średnie ciepło właściwe wody 0.00418 [GJ t-1 K-1].

 

O ile nie zostało to zaznaczone, dużymi literami alfabetu łacińskiego oznaczono wielomiany zmiennej z-1 np. A(z-1) = a0+a1z-1+...+anAz-nA, nA oznacza stopień wielomianu A.

 

B1+ - wielomian którego zera leżą w okręgu jednostkowym,

B1 - wielomian którego zera leżą na okręgu jednostkowym lub poza nim,

k,k1,k2,k3 – czasy opóźnień obiektu,

Ts – okres próbkowania.

 

 

1. Opis systemu ogrzewania w budynku B1

 

 

            Uproszczony schemat instalacji CO w budynku B1 przedstawiono na rys. 1. Instalacja ta zawiera dwa obwody. Obwód pierwszy (oznaczenia z indeksem 1) zasila halę technologiczną a drugi, zasila budynek B1 (oznaczenia z indeksem 2).

 

 

Rys.1. Uproszczony schemat instalacji CO w budynku B1.

 

 

Aktualne temperatury, przepływy i sterowania zaworami są mierzone co Ts = 10min i udostępnione po adresem http://wg.agh.edu.pl/b1.htm , historię pracy instalacji można zobaczyć na stronie http://wg.ia.agh.edu.pl/.

 

Każdy obwód składa się z kilkudziesięciu kaloryferów pracujących w różnych warunkach termicznych. Pojedynczy kaloryfer można zamodelować jako rurę z wodą umieszczoną w pomieszczeniu o temperaturze Tpom (rys.2.). 

 

 

Rys. 2. Uproszczony model kaloryfera.

 

 

Tego typu model przy pominięciu pojemności cieplnej ścian grzejnika, można opisać następującym równaniem różniczkowym cząstkowym

 

                                                  (1)

 

Współczynniki modelu l, A, B, są funkcjami przepływu, a współczynnik B zależy również od warunków wymiany ciepła w pomieszczeniu. Dodatkowo w modelu (1) należałoby uwzględnić moc wypromieniowywaną przez grzejnik. W stanie ustalonym temperatura wody powrotnej Tpco jest nieliniową rosnącą funkcją przepływu Fzco i zawsze zachodzi warunek

. Ponadto system (1) ma zmienne opóźnienie zależne od przepływu Fzco. Inne własności modelu (1) zob. [xx].

 

            Jakkolwiek zamodelowanie pojedynczego grzejnika nie stanowi obecnie poważnego problemu, to stworzenie bazującego na równaniach fizycznych modelu matematycznego  kilkudziesięciu kaloryferów i rur sprawia poważne trudności. Mamy bowiem do czynienia z systemem nieliniowym, niestacjonarnym, o kilkudziesięciu zmiennych opóźnieniach. Aby tego dokonać wykorzystamy modele typu ARMAX.

 

 

2. Modele instalacji CO

 

2.1 Charakterystyki zaworów

 

            Na rysunkach 3 i 4 przedstawiono zależność przepływu od stopnia otwarcia zaworów regulacyjnych z1 i z2. Charakterystykę zaworów przybliżono funkcją

 

.                                                                                                                   (2)

u – stopień otwarcia zaworu [%].

 

Parametr

Zawór  z1

Zawór  z2

A1

0.3939

0.58161

A2

10.22443

10.92176

u0

74.0276

50.08153

p

8.92594

4.10168

                                              

Rys. 3.  Zależność przepływu od stopnia otwarcia zaworu dla zaworu z1.

Rys. 4.  Zależność przepływu od stopnia otwarcia zaworu dla zaworu z2.

 

 

 

2.2 Model ARMAX instalacji CO

 

 

            Na podstawie analizy równań (1) można napisać następujący model dyskretny

 

 

                   (3)

e(i) - biały szum, i=1,2,... .

 

 

Celem identyfikacji była minimalizacja wskaźnika jakości

 

                                                                                 (4)

 

przy czym N oznacza liczbę próbek,  oznacza predykcję wyjścia obiektu wykonaną za pomocą równania (3) w chwili i–1 na chwilę i przy założeniu e(i)=0.  jest zmierzoną wartością temperatury powrotnej w chwili i. Nieznane wartości białego szumu estymujemy z zależności

 

.                                                                                             (5)

 

Powyższa metoda identyfikacji znana jest jako metoda predykcji błędu (PEM zob.[xx]), a odpowiedni algorytm identyfikacji jest zaimplementowany w przyborniku identification toolbox pakietu MATLAB. Identyfikacja została wykonana dla obwodu zasilającego budynek B1 (rys.1). Na rys. 6 przedstawiono odpowiedź obiektu i modelu oraz sygnały wejściowe obiektu.  Wyboru struktury obiektu dokonano kierując się minimalną wartością wskaźnika

 

                                                                                                          (6)

 

gdzie Np. oznacza liczbę parametrów modelu (3). Zauważmy, że w zależność (6)  nakłada karę za wybór modelu o zbyt złożonej strukturze (wybór struktury modelu zob. [xx s123]). Dla zadanych wartości nA , nB1 , nB2 , nB3 , nC, k1 , k2, k3, wykonywana była identyfikacja parametrów modelu (3), a następnie obliczano wartość wskaźnika (6), najlepiej dopasowany model miał strukturę

 

nA = 2, nB1 = 2, nB2 = 2, nB3 = 2, nC = 3, k1 = 2, k2 = 1, k3 = 6.

 

 

Równanie (3) można zapisać w postaci transmitancyjnej

 

 

                             (7)

przy czym

 

A = 1-1.731 z-1 + 0.753 z-2,

B1 = 0.2111 - 0.3024z-1 + 0.1 z-2,

B2 = 0.1171 - 0.1265z-1 + 0.0217 z-2,

B3 = 0.0564 - 0.0880 z-1 + 0.0392 z-2,

C = 1 - 0.3005 z-1  - 0.0827 z-2  - 0.0358 z-3,

 

Wielomiany A, B1, B2, B3, C są stabilne tzn. mają pierwiastki wewnątrz okręgu jednostkowego.

Rys. 6. Wynik identyfikacji dla modelu (3). Na wykresie nałożono temperaturę zasilania Tzco, przepływ Fzco, średnią temperaturę w pomieszczeniach Tpom, temperaturę Tpco oraz wynik jednokrokowej predykcji temperatury Tpco wykonanej w chwili i-1 na chwilę i. Błędy predykcji są lepiej widoczne na rysunku 7.

 

W celu weryfikacji jakości modelu (3) wykonano testy statystyczne dla reszt modelu. Na rysunku 8. Przedstawiono funkcję autokorelacji reszt modelu. Wykonano trzy różne testy statystyczne dla reszt modelu (3)(zob. [xx])

 

-         Test autokorelacji,

-         Test chi2,

-         Test znaku.

 

Poziom istotności wynosił 5%. Otrzymano następujące wyniki

 

 

t obliczone

t tablicowe

test autokorelacji

5.79

1.96

test chi2

98

55.7585

test znaku

3.33

1.96

 

 

Rys. 7. Powiększony fragment rysunku 6. Kolor czerwony – model, kolor czarny obiekt.

 

Rys. 8. Funkcja autokorelacji reszt modelu (3). Poziome linie odpowiadają granicom 95% przedziału ufności dla hipotezy, że reszty są białym szumem.

 

Jak widać reszty nie są białym szumem, jednak wartości obliczone i wartości tablicowe nie różnią się wiele, przyjmujemy zatem, że model (3) jest dobrze dopasowany. Odpowiedzi skokowe modelu (3) przedstawiono na rys.9

 

Rys. 9. Odpowiedzi skokowe modelu. Idąc od góry mamy odpowiedź na jednostkowy skok przepływu (Fzco), temperatury zasilania (Tzco), temperatury w pomieszczeniu(Tpom).

 

 

3. Model poboru mocy przez budynek B1

 

            Moc pobierana przez budynek w chwili i wynosi:

 

                                                                                               (8)

 

Na podstawie równania (3) i wzoru (8) wyznaczono jednokrokową predykcję mocy pobieranej przez budynek. Wyniki przedstawiono na rys. 10, 11.

 

Rys. 10. U góry moc pobierana przez budynek (kolor czarny) i jej predykcja jednokrokowa (kolor czerwony) wykonana w oparciu o model (3) i wzór (8). U dołu temperatury zasilająca, powrotna, w pomieszczeniu oraz przepływ. Błędy predykcji widoczne są lepiej na rysunku 11.

 

Rys. 11  Powiększony fragment rysunku 10. Kolor czarny moc pobierana, kolor czerwony predykcja jednokrokowa mocy pobieranej przez budynek.

 

Mamy zatem następujący model nieliniowy poboru mocy przez budynek

 

,                                                                                                         (9)

 

,                  (10)

 

.                                                                                  (11)

 

Poniżej przedstawiono realizację modelu w SIMULINKU.

 

 

 

Rys. 12. Model symulacyjny instalacji CO zrealizowany w SIMULINKU.

Kliknij na rysunek aby pobrać model i dane do identyfikacji.

 

 

 

Korzystając z (9), (10), (11) można wyznaczyć zależność mocy pobieranej przez budynek w stanie ustalonym w funkcji sterowania zaworem i temperatury wody zasilającej. Rys. 13, 14 pokazują tę zależność

 

Rys. 13. Zależność mocy pobieranej przez budynek B1 od stopnia otwarcia zaworu i temperatury zasilania. Przyjęto, że temperatura w pomieszczeniach jest stała i wynosi 20 °C.

 

Rys. 14. Zależność mocy pobieranej przez budynek B1 od stopnia otwarcia zaworu dla różnych  temperatur zasilania. Przyjęto, że temperatura w pomieszczeniach jest stała i wynosi 20 °C.

 

 

Uwaga 1

            Warto zauważyć pewną charakterystyczną cechę poboru mocy przez budynek. Skokowa zmiana przepływu Fzco powoduje natychmiastowy  wzrost mocy, która później powoli spada. Zjawisko takie występuje, gdyż temperatura wody powrotnej nie ulega zmianie natychmiast lecz rośnie dopiero po pewnym czasie co powoduje spadek mocy. Efekt ten znajduje potwierdzenie zarówno w danych pomiarowych jak i w symulacjach. Zob. rys. 10 i 15.

 

Rys. 15. Wynik eksperymentu. U góry moc dostarczana do budynku, u dołu przepływ Fzco. Skokowy wzrost przepływu powoduje szybki wzrost  mocy poczym następuje jej spadek.

 

 

 

3.1 Linearyzacja modelu poboru mocy przez budynek

 

            Niech Tpcou, Tzcou, Tpomu, Fzcou, uzcou, Pu oznaczają ustalone wartości wejść i wyjść modelu, oznaczmy

 

                                                                                                                                     (12)

 

Rozwijając (9) i (11) w szereg Taylore’a z dokładnością do wyrazów pierwszego rzędu

oraz korzystając z (7) otrzymamy

 

                 (13)

przy czym

 

                                                                                                             (14)

 

Odpowiedzi skokowe modelu (13) dla uzcou = 50%, Fzcou = 5.73 [t/h], Tzcou = 44°C, Tpomu = 20°C, Tpcou = 34.14°C

Przedstawiono na rysunku 16.

 

Rys. 16. Odpowiedzi skokowe modelu poboru mocy przez budynek . Idąc od góry mamy zmianę mocy w odpowiedzi na jednostkowy skok położenia zaworu (uzco), temperatury zasilania (Tzco), temperatury w pomieszczeniu(Tpom). Ustalone wartości wejść wynosiły uzcou = 50%, Fzcou = 5.73 [t/h], Tzcou = 44°C, Tpomu = 20°C, Tpcou = 34.14°C.

 

Zauważmy, że wzrost przepływu lub temperatury zasilania powoduje początkowy szybki skok mocy a następnie jej powolny spadek. Wzrost temperatury w pomieszczeniu powoduje zmniejszenie mocy dostarczanej do budynku.  Współczynniki wielomianów w liczniku transmitancji we wzorze (13) silnie zależą od ustalonych wartości wejść i wyjść modelu. Na Rys. 17,18 przedstawiono zależność współczynnika wzmocnienia od temperatury zasilania i położenia zaworu.

 

 

           

Rys. 17. Zależność współczynnika wzmocnienia od temperatury zasilania i położenia zaworu. Przyjęto, że temperatura w pomieszczeniach jest stała i wynosi Tpomu = 20°C.

Rys. 18. Zależność współczynnika wzmocnienia od położenia zaworu dla różnych temperatur zasilania. Przyjęto, że temperatura w pomieszczeniach jest stała i wynosi Tpomu = 20°C.

 

 

 

 

Tekst, rysunki, symulacje i eksperymenty: Piotr Bania,

Współpraca: Wojciech Grega, Krzysztof Kołek.