Článek přibližuje důvody a způsob použití moderních nástrojů pro matematické a fyzikální modelování Matlab® a Comsol Multiphysics® k nepřímé identifikaci tepelného výměníku jako významné komponenty systémů centrálního zásobování teplem. Představuje verifikovaný fyzikální model, jehož existenci nepřímá metoda identifikace předpokládá.
V současnosti je v popředí zájmu inteligentní dodávka energií, tj. tepla a elektřiny, spočívající v optimálním řízení okamžitého výkonu zdrojů s kombinovanou výrobou tepla a elektřiny a v řízení výroby, rozvodu a spotřeby energií inteligentními řídicími systémy. Teplo se přitom rozvádí prostřednictvím systémů centrálního zásobování teplem (SCZT), tj. geograficky rozlehlými sítěmi teplovodů a výměníkových stanic zajišťujícími dopravu tepla od jeho zdrojů na místa koncové spotřeby. Jde o systémy dynamické, u nichž se v čase mění skladba zdrojů i míst koncové spotřeby. Není pochyb o tom, že z hlediska technického, ekonomického i ekologického je celospolečensky žádoucí co možná nejvyšší kvalita řízení systémů CZT.
Inteligentní řízení SCZT s použitím modelů
Hospodárné řízení subsystému výroby a využití tepelné energie spočívá jednak v hospodárném rozdělování zatížení mezi jednotlivé spolupracující zdroje a mezi jednotlivé výrobní jednotky uvnitř těchto zdrojů a jednak v určení vhodné sestavy spolupracujících zdrojů, včetně určení vhodné skladby spolupracujících jednotek uvnitř zdroje, jakož i v určení vhodného ekonomicky podloženého spouštění či odstavování jednotlivých výrobních jednotek, popř. i celých zdrojů.
Ke zvládnutí uvedené úlohy je třeba vytvořit komplexní matematicko-fyzikální model rozlehlých teplárenských horkovodních i parních sítí a výrobních zdrojů tepla, který umožní najít nové nekonvenční řídicí algoritmy potřebné k optimálnímu – ze zmíněných hledisek technického, ekonomického a ekologického – řízení komplexního dodavatelského řetězce výroba, doprava + rozvod a spotřeba tepla zejména rozlehlých teplárenských soustav. Postup tvorby a způsob využití komplexního modelu SCZT jsou schematicky naznačeny na obr. 1. Cesta ke komplexnímu modelu vede přes modely jednotlivých prvků SCZT.
Systém CZT obsahuje, kromě zdroje tepla, celkem tři tyto základní skupiny prvků:
-
předávací výměníkové stanice,
-
primární a sekundární potrubní síť,
-
odběratelský segment, tj. systém topných těles (radiátorů).
Vlastnosti zmíněných základních prvků mají vliv na procesy distribuce a spotřeby tepelné energie. Stav jednotlivých prvků SCZT je třeba ovlivňovat (řídit) tak, aby bylo dosaženo požadovaného optimálního chování celého systému. K tomu je nutné znát stavy (hodnoty parametrů) jednotlivých prvků SCZT, resp. jejich odhady – což znamená prvky identifikovat.
Identifikovat je třeba především nejvýznamnější prvky (komponenty) SCZT z již uvedených základní skupin, tedy:
-
tepelný výměník jako hlavní komponentu výměníkové stanice,
-
potrubí,
-
typické topné těleso (radiátor ústředního topení).
Při identifikaci reálných komponent je třeba měřit, a v tomto ohledu se v rozlehlých SCZT v praxi naráží na značné potíže.
Zdroji těchto potíží jsou některé obecné vlastnosti (měřicích) systémů. Jde především o tyto jevy:
-
dynamické vlastnosti měřicího systému,
-
interakci měřicího systému s okolím (mechanická, tepelná a elektromagnetická kompatibilita),
-
nelinearitu měřicího systému,
-
historickou závislost, tj. změny vlastností měřicího systému v čase.
Východisko se zde nabízí v podobě použití nepřímé metody identifikace prostřednictvím fyzikálního modelu, který je verifikován s použitím modelu matematického.
Další text bude omezen na ukázku způsobu tvorby obou modelů v případě tepelného výměníku jako nejdůležitější komponenty SCZT z hlediska řízení systému. Fyzikální model výměníku byl vytvořen při použití programu Comsol Multiphysics a matematický (simulační) model při použití nástroje Matlab Identification Toolbox.
Poznamenejme, že při identifikaci a návrhu algoritmů řízení SCZT jako celku se z nástrojů obsažených v prostředí Matlab od firmy The MathWorks vedle již uvedeného nástroje Identification Toolbox dále uplatní nástroje Simulink a Control Toolbox, jak také ukazuje obr. 1.
Fyzikální model výměníku tepla
Fyzikální model je zdrojem údajů pro stanovení hodnot parametrů matematického modelu vybrané komponenty SCZT, v daném případě jednoduchého trubkového výměníku s geometrií podle obr. 2 a tab. 1.
K určení fyzikálního modelu výměníku tepla je potřebné nalézt řešení soustavy rovnic, které popisují proudění teplonosného média a přenos tepla uvnitř tohoto média v prostoru a v čase. Fyzikální děje ve výměníku tepla lze popsat s použitím Fourierovy-Kirchhoffovy rovnice a Navierovy-Stokesovy rovnice.
Fourierova-Kirchhoffova rovnice pro výměník tepla má tvar
rovnice (1)
kde
Cpje měrná tepelná kapacita látky při konstantním tlaku,
Q reprezentace vnitřních zdrojů: Q(x, y, z) = 0,
T reprezentace teplotního pole (skalárního),
k tepelná vodivost látky (obecně tenzor, jinak u homogenních látek skalár),
v reprezentace rychlostního (vektorového) pole,
t čas,
x, y, z prostorové souřadnice,
ρ hustota látky (např. teplonosné médium nebo materiál stěny výměníku).
Předpokládá se, že vnější plášť výměníku je ideálně izolovaný od okolí, takže platí skalární vztah
rovnice (2)
kde n0 je jednotkový normálový vektor.
Stěna výměníku oddělující primární okruh (médium 1) od sekundárního okruhu (médium 2) vyhovuje podmínkám kontinuity vedení tepla. Tato podmínka odpovídá rovnici
rovnice (3)
Okrajové podmínky na vstupu primárního a sekundárního média jsou
rovnice (4)
Okrajové podmínky na výstupu primárního a sekundárního média odpovídají okamžité hodnotě tepelného toku v daném místě, což odpovídá rovnici
rovnice (5)
Při přenosu tepla mezi primárním a sekundárním médiem v tepelném výměníku se dominantně uplatňuje nucená konvekce. K Fourierově-Kirchhoffově rovnici (1) je tudíž nutné připojit rovnici spojitosti a rovnice popisující proudění teplonosného média. V případě proudových polí s malými rychlostmi proudění je s ohledem na kapacitu paměti výpočetních prostředků přípustným kompromisem použít Navierovu-Stokesovu rovnici
rovnice (6)
kde
v je vektor rychlosti pohybu teplonosného média ve výměníku,
vx, vy, vzsložky vektoru rychlosti v,
F(x, y, z) síla působící na jednotku objemu tekutiny (F = 0).
Průběh vektoru rychlosti teplonosného média je počítán při zachování podmínky kontinuity
rovnice (7)
Předpokládá se, že vnější plášť i primární a sekundární potrubí výměníku jsou dokonale těsné a že teplonosné médium nikde neuniká, tj. jeho hmotnostní tok mezi primárním a sekundárním okruhem je konstantní a rovný nule. Z toho plyne okrajová podmínka, podle které rychlost pohybu média v této části výměníku je rovna nule, tedy
v (x, y, z) = 0 (8)
Okrajové podmínky na vstupu primárního a sekundárního média jsou určeny konstantními rychlostmi proudění tekutiny
v(t, x, y, z) = v11, v11 = 1 m·s–1
v(t, x, y, z) = v12, v12 = 2 m·s–1 (9)
Okrajové podmínky na výstupu primárního a sekundárního média odpovídají podmínkám při konstantním tlaku p média při nulovém viskózním napětí, tedy platí
p(t, x, y, z) = p21, p21 = 2,1 MPa
p(t, x, y, z) = p22, p22 = 1 MPa (10)
Řešení soustavy rovnic (1) a (6) při specifikovaných okrajových podmínkách (4) až (5) a (7) až (9) pro rozložení teploty v čase a prostoru pro výměník podle obr. 1 a tab. 1 je znázorněno na obr. 3. Rozložení rychlostního pole ukazuje obr. 4.
Simulační model výměníku tepla
Simulační model je založen na řešení matematického modelu v podobě stavové maticové lineární diferenciální rovnice ve tvaru
ds/dt = A(r) s(t) + B(r)ur(t) (11)
kde
s(t) je stavový vektor matematického modelu,
ur(t) vstupní veličina v podobě teploty na vstupu sekundárního okruhu,
r pracovní bod daný teplotou v primárním okruhu výměníku.
Pro výstupní veličinu matematického modelu, kterou je teplota média na výstupu sekundárního okruhu yT(t), platí
yT(t) = C(r)s(t) (12)
Matice A(r), B(r) a C(r) jsou závislé na parametru, který je určen pracovním bodem výměníku tepla. Pracovní bod výměníku je určen rychlostí toku v(t) a teplotou T(t) média v primárním okruhu. Matice matematického modelu (11) byly identifikovány s fyzikálním modelem, jehož pohyb byl vypočítán na základě řešení soustavy parciálních diferenciálních rovnic (1) a (6) s uvážením příslušných počátečních a okrajových podmínek (viz výše). Identifikace matematického modelu byla prováděna na dané množině pracovních bodů výměníku tepla. Výsledky identifikace jsou ukázány na obr. 5 a obr. 6, kde jsou zobrazeny průběhy výstupních teplot média ze sekundárního okruhu fyzikálního a matematického modelu a průběhy chyb identifikace v podmínkách hraničních pracovních bodů, kdy rychlost proudění média v primárním okruhu je v11(t) = 1 m·s–1 a teplota média v primárním okruhu je T11(t) = 353,15 K, popř. 390 K.
Matice matematického modelu určeného rovnicemi (11) a (12) pro pracovní bod r1 odpovídající hodnotám fyzikálních veličin primárního okruhu T11 = 390 K, v11 = 1 m·s–1 mají podobu
rovnice (13)
Obdobné jsou matice pro pracovní bod r2 odpovídající hodnotám fyzikálních veličin primárního okruhu T11 = 353,15 K, v11 = 1 m·s–1, které pro omezený rozsah článku již nejsou uvedeny.
Závěr
Funkce výměníku tepla je založena na pohybu tekutiny, který má nelineární charakter. Pro vypracování návrhu algoritmu řízení, popř. identifikace stavu výměníku tepla je potřebné navrhnout model popisující dynamické vlastnosti výměníku. Většina algoritmů řízení, popř. identifikace algoritmů využitelných v praxi je založena na lineární teorii řízení, popř. identifikace. Je proto vhodné nalézt takový model, který bude lineární v určité oblasti blízko pracovního bodu. Pohyb nelineárního fyzikálního modelu tak může být aproximován matematickým modelem, jehož vlastnosti jsou závislé na volbě parametru. Parametr matematického modelu je určen pracovním bodem fyzikálního modelu výměníku tepla. Předpokladem takové aproximace je verifikovaný fyzikální model představený v článku.
Poděkování
V článku jsou uvedeny výsledky získané v rámci řešení projektu MŠMT ČR Národního programu výzkumu II s názvem Informační technologie pro znalostní společnost (2C), evid. č. 2C06007.
Literatura:
[1] BROGAN, W. L.: Modern Control Theory. Prentice Hall Inc., 1991.
[2] REDDY, J. N. – GARTLING, D. K.: The Finite Element Method in Heat Transfer and Fluid Dynamics. CRC Press, 2000.
[3] HUEBNER, K. H. – DEWHIRST, D. L. – SMITH, D. E. – BYROM, T. G.: The Finite Element Method for Engineers. John Wiley & Sons Inc., 2001.
Ing. Jiří Marek, CSc.,
Obr. 1. Struktura procesu identifikace a řízení soustavy CZT
Obr. 2. Geometrie uvažovaného trubkového průtokového výměníku
Obr. 3. Rozložení teplot podle fyzikálního modelu (Comsol Multiphysics)
Obr. 4. Rozložení rychlosti podle fyzikálního modelu (Comsol Multiphysics)
Obr. 5. Časový průběh teploty na výstupu sekundárního okruhu fyzikálního T22F matematického modelu T22M absolutní chyby identifikace e(t) při pracovním bodu T11= 1 m·s–1= 390 K, v11(t) a(t) a
Obr. 6. Časový průběh teploty na výstupu sekundárního okruhu fyzikálního T22F matematického modelu T22M absolutní chyby identifikace e(t) při pracovním bodu T11= 1 m·s–1= 353,15 K, v11(t) a(t) a
Tab. 1. Rozměry uvažovaného trubkového výměníku