MATEMATICKÉ MODELY PŘENOSOVÝCH DĚJŮ PŘI PRODMÝCHÁVÁNÍ OCELI V LICÍ PÁNVI INERTNÍM PLYNEM

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________

MATEMATICKÉ MODELY PŘENOSOVÝCH DĚJŮ PŘI PRODMÝCHÁVÁNÍ OCELI V LICÍ PÁNVI INERTNÍM PLYNEM MATHEMATICAL MODELS OF TRANSITION PHENOMENA AT GAS BLOWING INTO LADLE Jan Morávkaa Karel Michalekb a

Třinecký inženýring, a.s., Frýdecká 126, 739 61 Třinec - Staré Město, ČR, [email protected] b VŠB - Technická univerzita Ostrava, 17.listopadu 15, 708 33 Ostrava - Poruba, ČR, [email protected] Abstrakt Teplotní, chemická a fyzikální homogenizace oceli v licí pánvi se v praxi uskutečňuje prostřednictvím dmýchání inertního plynu (argonu) do této lázně. Příspěvek popisuje fyzikální a kybernetický přístup k hledání vhodného aproximačně-regresního modelu naměřené normované koncentrace stopovací látky ve zmenšeném fyzikálním modelu licí pánve při prodmýchávání oceli inertním plynem. Procesy probíhající při dmýchání plynu lze aproximovat pomocí Heavisideova jednotkového skoku a odezva koncentrace má pak charakter přechodové funkce. Fyzikální přístup umožňuje sestavit adekvátní matematický model děje ve tvaru tzv. bílé skřínky, kde je známá struktura (i parametry) modelu. Kybernetický přístup vychází pouze z naměřených vstupů a výstupů (tzv. černá skřínka), případně i z doplňkových podmínek (tzv. šedá skřínka) a strukturu modelu podle nich volí. Tyto modely se označují jako empirické modely. V příspěvku jsou prezentovány a srovnány čtyři modely – jeden fyzikálně adekvátní a tři empirické. Abstract Temperature, chemical and physical homogenisation of steel in the casting ladle is practically done by blowing inert gas (argon) into the bath. The paper describes the physical and cybernetic approach to the task of looking for an appropriate approximation-regression model of the measured standardised concentration of the tracer at a scale physical model of the ladle when the bath is blown by inert gas. The processes occurring while the gas is being blasted may be approximated using the Heaviside unit step and the reaction of the concentration then has the character of a transition function. The physical approach allows assembling an adequate mathematical model of the processes in the shape of a so-called white box, where the structure of the model are known. The cybernetic approach only draws on the measured inputs and outputs (so-called black box) as well as any additional conditions (socalled grey box) and the structure of the model is chosen according to them. These models are referred to as empirical models. The paper presents and compares four models – physically adequate one and three empirical ones. 1. ÚVOD Pro naměřené časové průběhy koncentrací stopovací látky ve fyzikálním modelu licí pánve (LP) při prodmýchávání oceli inertním plynem (argonem) bylo vhodné sestavit fyzikálně matematické (fyzikálně adekvátní) a empirické modely. Čtyři navržené modely byly verifikovány pomocí nelineární regresní analýzy a jejich výstupy byly porovnány.

1

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________ Výsledky analýz mohou být použity pro nastavení vhodného režimu práce prodmýchávání, jako i pro výuku na technických univerzitách. 2. POPIS SITUACE Schematické znázornění situace při prodmýchávání oceli (vody) inertním plynem (argonem) v modelu licí pánve (označený dále mLP, který byl vytvořený ve zmenšeném geometrickém měřítku 1:10) je na mLP obr.1: hk

K1

Obr. 1. Dmýchání argonu do modelu pánve Fig. 1. Blowing of argon into the ladle model

I q

V nádobě modelu LP jsou zasunuta tři vodivostní čidla K1, K2 a K3. Z excentricky II (mimosově) umístěného prodmýchávacího prvku P (dmýchací tvárnice, tzv. „kamen“) na K3 P dně mLP proudí bubliny argonu o konstantním objemovém průtoku q. Tyto „rozrušují“ vrstvu koncentračně obohacené a zabarvené vody o tloušťce hk a dochází k postupnému promíchávání obohacené a čisté kapaliny (oceli, vody). Na molekuly vody v blízkosti čidel v podstatě působí proti sobě dvě (tlakové) síly I a II.

K2

h

3. NAMĚŘENÁ DATA Analýza a syntéza fyzikálně matematického modelu byla uskutečněna na základě experimentu ze dne 12.5.2006, kde průběhy naměřených koncentrací ci(t), i ∈ {1,2,3}, (s periodou snímání ∆t ≈ 0.5 s), na čidlech K1, K2 a K3 měly následující tvar - obr.2: 0.09

K1

K2

K3

c(t) - koncentrace [% hm.]

0.08

K1

0.07 0.06 0.05 0.04 0.03

cu 0.02

K2 K3

cp

0.01 0.00 0

10

20

30

40

50

60

70

80

90

100 110 120 130 140 150 160 170 180 190 200

t - čas [s]

Obr. 2. Průběhy koncentrací na čidlech K1, K2 a K3 Fig. 2. Development of concentration in the sensors K1, K2 and K3

•

Z průběhů koncentrací na čidlech je zřejmých několik skutečností: zahájení a průběh prodmýchávání oceli plynem lze aproximačně uvažovat ve tvaru Heavisideova jednotkového skoku a tím lze průběhy koncentrací chápat jako přechodové

2

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________ charakteristiky, • čidla zareagovala až po uplynutí určitého dopravního zpoždění, které je úměrné vzdálenosti jednotlivých čidel od hladiny lázně v mLP, • překmit průběhů (zřejmě úměrný velikosti síly I nebo spíše rozdílu sil I a II) také klesá se vzdáleností od hladiny lázně v mLP, • ustálená (konečná) hodnota koncentrace je úměrná podílu objemů koncentračně obohacené a čisté modelové kapaliny (vody). Jak je vidět z obr.2, počáteční koncentrace není nulová, ale odpovídá zbytkové (přirozené) vodivosti modelové kapaliny (vody). Z tohoto důvodu, ale i z důvodu porovnání průběhů na všech čidlech, je vhodné zavést normovanou (a bezrozměrovou) koncentraci podle vztahu: cn (t ) =

c(t ) − c p cu − c p

(1)

[−]

kde je

t - čas [s], c(t) - naměřená koncentrace zabarvené vody [hm. %], cp - počáteční hodnota koncentrace [hm. %], cu - ustálená (konečná) hodnota koncentrace [hm. %]. Ze vztahu (1) je zřejmé, že počáteční hodnota normované koncentrace bude nulová a konečná (ustálená) hodnota bude rovná jedné. V důsledku vzniku místních koncentračních heterogenit v lázni mohou hodnoty této normované koncentrace v průběhu dmýchání překročit hodnotu jedné. 4. FYZIKÁLNĚ MATEMATICKÝ MODEL Přesný teoretický popis dějů při prodmýchávaní oceli argonem by byl velice složitý a vedl by k soustavě nelineárních parciálních diferenciálních rovnic popisujících přenos hybnosti, tepla, složek a s budicí funkcí ve tvaru rovnice tzv. deterministického chaosu (probublávání argonu). Takovýto systém by byl řešitelný pouze numericky, např. prostřednictvím CFD (Computational Fluid Dynamics) programů (např. programu Fluent). Na základě schématu mLP a průběhů koncentrací však byl navržen zjednodušený, tzn. minimální, nejjednodušší možný (vycházející z použití principu tzv. Occamovy břitvy: „modely nemají být složitější než je nezbytně nutné“) lineární fyzikálně matematický, či jinak fyzikálně adekvátní model chování koncentrace oceli v LP při jejím prodmýchávání. Popisovaný děj byl chápán ve tvaru kybernetického modelu, který je možné přehledně znázornit pomocí tzv. blokového GI schématu na obr.3: Obr. 3. Blokové schéma mLP Fig. 3. Block diagram of the ladle model

Q(s)

Gd

+ GII

C(s)

-

(model platí pro libovolné čidlo K1, K2 nebo K3, označení Q(s) znamená Laplaceův (L-)obraz průtoku plynu, C(s) je L-obraz koncentrace). Jde tedy o sériovo-paralelní zapojení tří členů, a to členu dopravního zpoždění a dvou paralelně a proti sobě působících proporcionálních (setrvačných) členů. Přenos dopravního zpoždění má tvar: Gd (s) = exp(−Td ⋅ s) = e −Td ⋅s [−]

kde je

t c(t)

(2)

- čas [s], - naměřená koncentrace zabarvené vody [hm. %],

3

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________ cp cu

- počáteční hodnota koncentrace [hm. %], - ustálená (konečná) hodnota koncentrace [hm. %].

Pro část modelu bez dopravního zpoždění předpokládáme dvě paralelně a protichůdně zapojené nejjednodušší proporcionální soustavy se setrvačností 1. řádu (s přenosy GI a GII), které jsou používány i v (chemické) kinetice procesů (zde však v ekvivalentním tvaru obyčejných diferenciálních rovnic 1. řádu s konstantními koeficienty): k1 , T1s + 1 k2 , GII (s) = T2 s + 1

(3a)

GI (s) =

kde jsou

(3b)

k1, k2 - koeficienty přenosu (zesílení) soustav [%⋅s/m3], T1, T2 - časové konstanty soustav [s].

Na základě analogie s tzv. rozdílovým termočlánkem [1] a na základě časových průběhů naměřených koncentrací (obr.2) je možné předpokládat výrazný (řádový) rozdíl v hodnotách časových konstant, tj. T2 >> T1 (čili sestupná část přechodové charakteristiky má podstatně větší časovou konstantu než část náběhová) a obecně také nestejné hodnoty koeficientů přenosu k1 ≠ k2. Pro uvedenou část modelu je potom možné na základě algebry přenosů sestavit následující spojité L-přenosy (pro nulové počáteční podmínky): G(s) = GI ( s) − GII ( s) , (k − k ) + (k1T2 − k2T1 )s k1 k , G( s ) = − 2 = 1 2 T1s + 1 T2 s + 1 (T1s + 1)(T2 s + 1)

(4) (5)

odkud pro L-obraz a originál (získaný pomocí zpětné L-transformace) přechodové funkce H(s) a h(t) (odezva na Heavisideův jednotkový skok průtoku inertního plynu, Q(s) = 1/s), jako i pro její mezní hodnoty dostaneme: 1 H (s) = C (s) = G(s) ⋅ Q( s) = G(s) ⋅ , s h(t ) = (k1 − k2 ) − k1 ⋅ exp(−t / T1 ) + k2 ⋅ exp(−t / T2 ) , h(0) = (k1 − k2 ) − k1 + k2 = 0, h(+∞) = k1 − k2 .

(6) (7) (8)

Jelikož ustálená hodnota normované koncentrace je rovná jedné, platí pro hodnotu přenosového koeficientu k2 vztah (vyjadřující také relaci, že k2 < k1): h(+∞) = k1 − k2 = 1 ⇒ k2 = k1 − 1 .

(9)

Dosazením uvedené relace do vztahu (7) dostaneme konečný výraz pro přechodovou funkci normované koncentrace, který může být současně použit jako nelineární regresní model se třemi parametry k1, T1 a T2: h(t ) = cn (t ) = 1 − k1 ⋅ exp(−t / T1 ) + (k1 − 1) ⋅ exp(−t / T2 ) .

(F1)

Pro úplnost sestavme ještě L-přenos této části soustavy, kdy po dosazení získaných relací pro přenosový koeficient do vztahu (5), dostaneme: G( s ) =

{k (T2 − T1 ) + T1}s + 1 ,

(10a)

(T1s + 1)(T2 s + 1) TD s + 1 G( s ) = , TD = f (k , T1 , T2 ) , (T1s + 1)(T2 s + 1)

(10b)

4

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________

TD = k (T2 − T1 ) + T1 ,

(10c)

odkud je zřejmé, že jde o modifikaci tzv. reálného derivačního členu se zpožděním 2.řádu, přičemž derivační časová konstanta TD je funkcí všech tří parametrů. Tento přenos (a chování modelu) ještě lépe odpovídá přenosu reálného PD regulátoru, avšak s navzájem různými reálnými póly (klasický reálný PD regulátor má dva póly komplexně sdružené) – viz [2]. K uvedenému modelu F1 se lze dostat také jinou cestou, a to logickou úvahou. Křivky podobného typu lze poskládat z exponenciál či použít rovnice chemické kinetiky (což může být téměř či úplně totéž). Nejjednodušší je součet dvou exponenciál (rostoucí a klesající) a konstanty: c(t ) = a1{1 − exp(−a2t )} + a3 ⋅ exp(−a4t ) + a5 ,

(11)

který vede při podmínkách c(0) = 0, c(∞) = 1 a po formálním přiřazení b1 = a1, b2 = a2, b3 = a4, na regresní rovnici (F1m), odpovídající modelu F1, přičemž mezi jejich koeficienty platí ekvivalence b1 = k1, b2 = 1/T1 [1/s], b3 = 1/T2 [1/s]: cn (t ) = 1 − b1 ⋅ exp(−b2t ) + (b1 − 1) ⋅ exp(−b3t ) .

(F1m)

Koeficienty b2 a b3 mají přitom v uvedené rovnici charakter jakýchsi rychlostních (či frekvenčních) konstant směšovacího děje. 5. EMPIRICKÉ MATEMATICKÉ MODELY Kromě získaného fyzikálně adekvátního matematického (a jemu odpovídajícího regresního) modelu F1, či modifikovaného F1m, byly navrženy ještě tři empirické modely (nazývané také empirické vzorce, funkce, vztahy, či formule) – typy a principy jejich výběru jsou uvedeny např. v [3], [4], [5]. Rozdíl mezi oběma typy modelů spočívá ve skutečnosti, že fyzikálně adekvátní model (zvaný též model teoretický, deterministický, či fenomenologický) odpovídá (i když většinou zjednodušeně) fyzikálním zákonům, jeho parametry mají fyzikální smysl a lze jej tedy použít i pro extrapolaci naměřených hodnot. Empirické matematické (regresní) modely tyto vlastnosti obecně nemají a „snaží“ se pouze čím jak nejlépe postihnout průběhy (trendy) v datech, přičemž je správné jejich použití pouze pro interpolaci hodnot (uvnitř intervalu nezávislé proměnné). Empirické modely, splňující mezní podmínky cn(t = 0) = 0, cn(t → +∞) = 1 pro normované koncentrace, mají následující tvar: cn (t ) =

t b1 − b2 +1 , b2 + b3 ⋅ t b4

(E1)

cn (t ) =

b1t + b2t 2 + b3t 3 . 1 + b4t + b5t 2 + b3t 3

(E2)

Obě empirické funkce mají tvar racionálně lomených (mocninných) funkcí, přičemž model E1 má 4 parametry a model E2 má dokonce 5 parametrů. Další možností je empirický model sestavený jako součet dvou funkcí: Hoerlovy funkce (f1) a vhodného násobku funkce arkus tangens (f2), přičemž pro jejich mezní hodnoty platí: f1(0) = f2(0) = 0, f1(∞) = 0, f2(∞) = 1: cn (t ) = b1 ⋅ t b2 ⋅ exp(−b3t ) +

2

π

⋅ arctg(b4t ) .

(E3)

O uvedených empirických matematických modelech lze říci, že mají charakter (z

5

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________ kybernetického pohledu) tzv. šedé skřínky, tzn. kromě vstupu a výstupu byly známy ještě doplňkové (vedlejší) podmínky [5] ve tvaru výše uvedených mezních podmínek. 6. POROVNÁNÍ MODELŮ Pro data LP_K1n.dat normované koncentrace, měřené na čidle K1, byly porovnány výsledky, získané pro čtyři modely, pomocí nelineární regrese. Výsledky jsou uvedeny v tab.1 (kde kritérium SSE = Sum of Squared Error, součet čtverců odchylek): Tabulka 1. Výsledky nelineární regresní analýzy pro data LP_K1n Table 1. Results of non-linear regression for data LP_K1n k1 T1 T2 ISE R2 Model b4 b5 b1 b2 b3 SSE [%] F1 6.238 3.214 26.633 19.33 94.12 E1 0.252 0.416 0.000057 2.751 22.15 93.26 E2 3.058 -0.023 0.00062 0.578 -0.0041 18.42 94.40 E2r 2.791 -0.018 0.00050 0.490 18.48 94.38 E3 1.743 0.512 0.0558 0.605 17.92 94.55

Poznámka 3 par. 4 par. 5 par. 4 par. 4 par.

Grafické znázornění naměřené a následně normované koncentrace na čidle K1, jako i regresní průběhy čtyř modelů, jsou viditelné na obr.4: Porovnání modelů 6 cn1(t)

F1

5

F1 E1

cn1 [-]

4

E2 3

E3

2 1 0 0

20

40

60

t [s]

80

Obr. 4. Porovnání regresních modelů pro data LP_K1n.dat Fig. 4. Comparison of the regression models for the data LP_K1n.dat

• •

Ze porovnání je zřejmé, že: modely jsou kvalitativně srovnatelné (s ohledem na „integrální“ kvalitu aproximace, vyjádřenou ukazateli R2 a SSE), přičemž R2 je vysoké, a to v rozmezí 93.3 ÷ 94.6 %, u empirického modelu E2 došlo v rozmezí hodnot času asi 131.5 až 108 sekund k „podkmitu“ pod hodnotu 1.0 a koeficient b5 je statisticky nevýznamný, takže by stačil modifikovaný a redukovaný model se 4 parametry (E2r): cn (t ) =

b1t + b2t 2 + b3t 3 , 1 + b4t + b3t 3

(E2r)

6

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________

•

•

který již má všechny koeficienty statisticky významné a index determinace R2 u něj klesl pouze o 0.02 %, obdobně u empirického modelu E3, který má prakticky shodný průběh s modelem E2, došlo v rozmezí hodnot času asi 144 až 107 sekund k „podkmitu“ pod hodnotu 1 a koeficient b4 je statisticky nevýznamný (v tomto modelu jej však nelze vypustit), nicméně nejlepší je model F1, který je fyzikálně adekvátní, má nejmenší počet parametrů a v počátečním náběhu poskytuje přiměřený nárůst hodnot (empirický model E1 má zde nárůst hodnot nepřiměřeně strmý).

7. REGRESNÍ MODEL U VŠECH KONCENTRACÍ V tab.2 je viditelné použití modelu F1 pro normované hodnoty koncentrace na všech třech čidlech (bez uvažování dopravního zpoždění): Tabulka 2. Výsledky nelineární regresní analýzy pro model F1 Table 2. Results of non-linear regression for the model F1 Čidlo k1 T1 T2 R2 [%] Poznámka K1 6.2 3.2 26.6 94.12 K2 1.4 2.9 45.1 74.98 min. R2 K3 1.1 5.7 46.5 96.62 max. R2 Na obr.5 jsou zobrazeny průběhy regresní funkce, vycházející z modelu F1, pro normované hodnoty koncentrace na všech třech čidlech K1, K2 a K3: Výsledky regrese pro tři čidla K1, K2 a K3 1.5

6

1.25

4

1 K3

3

0.75

2

0.5

1

cn2,3 [-]

cn1 [-]

K2 5

0.25 K1

0

0 0

50

100

t [s]

150

Obr. 5. Průběhy regrese modelu F1 pro data ze všech čidel Fig. 5. Courses of regression of the model F1 for the data from all three sensors Z výsledků je jasné, že model F1 je vhodný a použitelný pro popis průběhu koncentrace u všech tří čidel. Nejhorší přiblížení, tj. s R2 nižším asi o 20 % vůči ostatním čidlům, poskytovala regrese údajů koncentrace z prostředního čidla K2. Je také vidět, že koeficient přenosu (zesílení) k1 klesá se vzdálenosti čidla od hladiny lázně v mLP, zatímco časová konstanta T2 u této závislosti stoupá. 8. POUŽITÍ MODELU PRO HODNOCENÍ EXPERIMENTŮ Je zřejmé, že navržený fyzikálně adekvátní model může být použit pro hodnocení výsledků fyzikálního modelování dějů při prodmýchávání oceli v modelu licí pánve. U provedené série

7

METAL 2007 22.-24.5.2007, Hradec nad Moravicí ___________________________________________________________________________ (asi 120) experimentů v modelové laboratoři fyzikálního modelování v Třineckých železárnách, a.s. byly měněny některé základní veličiny simulace, a to: • q – objemový průtok dmýchaného argonu, • e – excentricita, poloha umístění prvků prodmýchávání, • n – počet prvků prodmýchávání. Obecně může být měněn i: • V – objem lázně v modelu LP, • ρ – hustota „horní“ a „spodní“ vrstvy lázně (ρh, ρs) • c – a jejich koncentrace(ch, cs). Koncentrace byla v průběhu měření navíc snímaná na třech čidlech ci(t), i ∈ {1,2,3}, což je nutné také zohlednit. Z uvedeného vyplývá, že parametry, či „konstanty“ fyzikálně adekvátního modelu Td, k1, T1 a T2 jsou vlastně funkcemi vstupních parametrů simulace, tzn. platí: Td = Td(V, ρh, ρs, ch, cs, q, e, n, i) atd. Díky normování koncentrací (ch, cs), použití stejné lázně (ρh, ρs) o stejném objemu (V), jako i jednoho dmýchacího elementu (n = 1) v určité poloze (e) lze předpokládat, že parametry modelu by byly funkcemi pouze dvou parametrů simulace, tj. že by platilo: Td = Td(q, i). Vyhodnocení závislostí parametrů modelu na vstupních parametrech provedených simulací je obsahem dalšího příspěvku. 9. ZÁVĚR Pro naměřené koncentrace zabarvené lázně v modelu licí pánve byly navrženy a ověřeny čtyři matematické (regresní) modely, a to jeden fyzikálně adekvátní a tři empirické. Jako nejvhodnější a přitom nejjednodušší (co odpovídá principu tzv. Occamovy břitvy: „nejjednodušší je obvykle správné, či alespoň vhodné“) z uvažovaných se ukázal být fyzikálně adekvátní model, který může být (a také byl) použit pro další analýzu vlivu vstupních parametrů simulace na koeficienty tohoto modelu, jako i pro nastavení vhodného (optimálního) režimu práce prodmýchávání oceli v licí pánvi. Nezanedbatelná není ani didaktická použitelnost obsahu příspěvku pro výuku studentů na (technických) univerzitách. Práce vznikla v rámci řešení grantového projektu č.106/07/0407 za finanční podpory Grantové agentury České republiky. LITERATURA [1] VÍTEČEK, A. & SMUTNÝ, L. & KUSYN, J. 1988. Teorie řízení I. Základní pojmy a řešené příklady z Laplaceovy transformace. III. vyd. Ostrava : skripta FSE VŠB, 1988. 52 s. [2] KUBÍK, S. aj. 1974. Teorie regulace - I. Lineární regulace. Praha : SNTL, 1974. 272 s. [3] BRONŠTEJN, I.N. & SEMENĎAJEV, K.A. 1964. Príručka matematiky pre inžinierov a pre študujúcich na vysokých školách technických. III. vyd. Bratislava : SVTL, 1964. 748 s. [4] PECHOČ, V. 1981. Vyhodnocování měření a početní metody v chemickém inženýrství. 2.přepr. vyd. Praha : SNTL, 1981. 228 s. [5] KROPÁČ, O. 1987. Náhodné jevy v mechanických soustavách. Praha : SNTL, 1987. 408 s.

8