Matematikai m

1
Matematikai módszerek a reakciókinetikában

• Turányi Tamás
• ELTE Kémiai Intézet
• Reakciókinetikai Laboratórium

2006
2
Matematikai módszerek a reakciókinetikában
Mérési adatok feldolgozása (alkalmazott
statisztika) ? Keszei, Szepesváry
Numerikus módszerek reakció-diffúzió rendszerek
szimulációjára ? Lagzi-Izsák
Részletes összetett reakciómechanizmusok
vizsgálata Reakcióutak vizsgálata Idoskála-ana
lízis Érzékenység-analízis Reakciómechanizmuso
k redukciója Bizonytalanság-analízis
3
1. Reakciókinetikai alapismeretek
4
1 Reakciókinetikai alapismeretek

• Kémiai változások leírása sztöchiometriai
  (bruttó) egyenlettel
• megmutatja a reaktánsok és termékek arányát
• általában nem játszódik le fizikailag

2 H2 O2 2 H2O 0 - 2 H2 -1
O2 2 H2O
?1 - 2 A1 H2 ?2 - 1 A2 O2 ?3
2 A3 H2O
sztöchiometriai együttható (reaktánsra negatív,
termékekre pozitív elemek) Az anyagfajták
sorrendje tetszoleges A sztöchiometriai
együtthatók szorzótényezo erejéig határozatlanok

5
Reakciósebesség
koncentrációváltozási sebesség
reakciósebesség
Yj az Aj moláris koncentrációja mól dm-3

kis koncentrációtartományban mindig igaz

k reakciósebességi együttható
?j reakciórend a j-edik anyagfajtára
bruttó reakciórend
6
Összetett reakciómechanizmusok
Szinte mindig sok reakciólépés játszódik le
egyszerre
A reakciólépések lehetnek elemi reakciók, amik
fizikailag így játszódnak le, de lehetnek több
elemi reakció összevonásából kapott reakciólépések
bal oldali sztöchiometriai együtthatók mátrixa
elemi reakció összegük legfeljebb 2 nulla vagy
pozitív egész szám nem elemi reakció nulla vagy
tetszoleges pozitív egész szám
jobb oldali sztöchiometriai együtthatók
mátrixa elemi reakció összegük legfeljebb 2
nulla vagy pozitív egész szám nem elemi reakció
tetszoleges valós szám (lehet nulla,
negatív, tört stb.)
sztöchiometriai együtthatók mátrixa
(INFORMÁCIÓVESZTÉS !!!)
7
Kinetikai differenciálegyenlet-rendszer
tömeghatás törvénye (Guldberg és Waage, 1865)

ki i-edik reakciólépés reakciósebességi
együtthatója
ri i-edik reakciólépés sebessége
Kinetikai differenciálegyenlet-rendszer
Kinetikai differenciálegyenlet-rendszer
mátrix-vektor alakban
8
Kinetikai differenciálegyenlet-rendszer egy
példa
Belouszov-Zsabotyinszkij oszcilláló reakció
(egyik) Oregonátor modellje 1. X Y ? 2
P k1 r1 k1xy 2. Y A ? X P k2 r2
k2ya 3. 2 X ? P A k3 r3 k3x2 4. X A ? 2
X 2 Z k4 r4 k4xa 5. X Z ? 0,5 X
A k5 r5 k5xz 6. Z M ? Z Y k6 r6 k6zm
X HBrO2 Y Br? Z Ce4 A BrO3? P HOBr M
malonsav
80 reakciólépést tartalmazó mechanizmust lehetett
redukálni erre a 6 reakcióra. Negatív és tört
jobb oldali sztöchiometriai együtthatók!
9
Kinetikai differenciálegyenlet-rendszer egy
példa 2.
1. X Y ? 2 P 2. Y A ? X P 3. 2 X ? P
A 4. X A ? 2 X 2 Z 5. X Z ? 0,5 X A 6. Z
M ? Z Y
X HBrO2 diffegyenlettel számítjuk Y
Br? diffegyenlettel számítjuk Z
Ce4 diffegyenlettel számítjuk A
BrO3? állandó koncentrációjú P HOBr csak
termék M malonsav állandó koncentrációjú
?
?
?
A kezdeti érték probléma
10
Jacobi-mátrix számítása



11
F-mátrix számítása



12
Jacobi-mátrix sajátérték-sajátvektor felbontása
Jacobi-mátrix nagyon fontos, mert megadja,
hogyan változik fi ha yj megváltozik! J csak
elsorendu fogyasztó és nulladrendu termelo
reakciókat tartalmazó rendszerekben állandó
(egyébként nem!).
Jacobi-mátrix sajátérték-sajátvektor felbontása
sajátértékeket tartalmazó diagonális mátrix
(komplex sajátértékek!)
baloldali sajátvektorok (sorvektorok) mátrixa
jobboldali sajátvektorok (oszlopvektorok) mátrixa
A bal- és jobboldali sajátvektorok ortonormáltak
Ennek következtében
13
Merev differenciálegyenlet-rendszerek
Jacobi-mátrix sajátérték-sajátvektor
analízise merevség Re(?1)/Re(?n) Nagyon
merev (stiff) diffegyenletrendszer ? csak
különleges algoritmussal oldható meg ? a
lépéshossz változik a merevséggel
14
Kinetikai diffegyenletrendszer tulajdonságai

• Egy anyag koncentrációváltozási sebességére
  felírt differenciálegyenlet-rendszer csak
  elsorendu deriváltat tartalmaz, ami a
  koncentrációk nem feltétlenül lineáris függvénye.
• ? elsorendu nemlineáris differenciálegyenlet-ren
  dszer
• Mivel általában minden anyag több reakcióban is
  részt vesz, ezért anyagok koncentrációváltozásai
  erosen csatoltak.
• ? erosen csatolt
• A reakciók sebessége igen sok nagyságrendet
  átfog. (égéskémiában 8, légkörkémiában 24)
• ? erosen merev
• Légkörkémiai modellekben a reakciók sebességi
  állandói explicit módon függhetnek az idotol,
  nyomás, a homérséklet és a sugárzási viszonyoktól
  való függésük miatt.
• ? lehet autonóm és nem autonóm
• A valós folyamatok többnyire térben inhomogének,
  így a kémiai folyamatok mellett transzport és más
  fizikai folyamatokkal is számolni kell.
• ? parciális diffegyenlet-rendszer, kémiai
  forrástaggal

15
Trajektória
Zárt kinetikai rendszer A koncentrációk addig
változnak, amíg a (termodinamikai) egyensúlyi
pontba nem érnek. Nyílt kinetikai
rendszer Folyamatosan adagoljuk a reaktánsokat
és eltávolítjuk a termékeket. Pl. jet-stirred
gas reactor, PSR perfectly stirred reactor CSTR
continuous stirred tank reactor, continuous-flow
stirred tank reactor Trajektória a megoldás
pályája a koncentrációtérben Zárt kinetikai
rendszer trajektóriája Kezdeti koncentrációk ?
egyensúlyi pont nyílt kinetikai rendszer
trajektóriája Kezdeti koncentrációk ?
stacionárius pont Kezdeti koncentrációk ? zárt
görbe (határciklus)
16
Megorzött tulajdonságokconserved properties
Izolált rendszer Az összes entalpia
állandó Zárt kinetikai rendszerben a
koncentrációk összege állandó, ha minden
reakciólépés mólszám-megorzo (igaz formális
rendszerekre is!) Zárt kinetikai rendszerben,
elemi reakciók esetén Az elemek száma
állandó Atomcsoportok (pl. benzolgyuru) száma
állandó lehet Megorzött tulajdonság C-atom ? 2
C2H4 1 CH4 6 C6H6 állandó H-atom ? 4
C2H4 4 CH4 6 C6H6 állandó A
koncentrációk lineáris konmbinációja
állandó Ahány megorzött tulajdonság, annyi ilyen
összefüggés!
N megorzott tulajdonság A sztöchiometriai mátrix
rangja N-el csökken ? N-el kevesebb (n-N)
változóval is pontosan szimulálható minden
koncentráció!
17
A reakciókinetika gyakorlati jelentosége

• Légkörkémiai folyamatok modellezése
• légszennyezés elorejelzése (idojárás elorejelzést
  igényel)
• kibocsátási korlátok megállapítása
• Égések, robbanások modellezése
• Pl. eromuvek, kazánok, motorok
• hatásfok optimalizálása
• szennyezoanyag-kibocsátás csökkentése
• Vegyi üzemek, gyártási folyamatok modellezése
• hatásfok és környezetvédelem optimalizálása
• Biokémiai folyamatok modellezése (systems
  biology)
• Metabolizmus hálózatok (pl. gyógyszerlebomlás
  leírása)
• Molekuláris jelterjedés
• Sejtciklus modellezése
• Reakciókinetikai formalizmussal leírt nem kémiai
  modellek

18
2. Reakcióutak
19
2. Reakcióutak
Reakcióutak milyen anyagból milyen anyag lesz?
Reakciófluxusok A nyilak vastagsága arányos az
átalakulás sebességével
Általában a könyvekben nem árulják el, hogyan
készültek a reakciófluxus ábrák. Warnatz J.,
Maas U., Dibble R. W. Combustion. Physical and
chemical fundamentals, modeling and simulation,
experiments, pullutant formation Springer, New
York, 1996 tüzeloanyagban gazdag
metán-levego láng
20
Reakció fluxusok
S.R. Turns An introduction to combustion.
Concepts and applications. second edition,
Boston, McGraw-Hill, 2000.
minden nyíl egy reakciólépést jelöl a nyíl
vastagsága arányos a reaktáns fogyási
sebességével
Nem jó ötlet, mert különbözo vastagú egymás utáni
nyilakhoz azonos fluxusok tartozhatnak ?
Megorzött változó fluxusát kell ábrázolni
!!! ? elemfluxusok (Revel et al., 1994)
21
Elemfluxusok
CH3 C3H7 gt C4H8 H2
reakciósebesség r1 H-atomok száma 3
7 8 2   H-atomok száma a bal oldalon
10   H-atomok fluxusa az egyik anyagfajtáról
a másikra CH3 ? C3H7 0 CH3 ? C4H8 3/108r1
2.4r1 CH3 ? H2 3/102r1
0.6r1 C3H7 ? CH3 0 C3H7 ? C4H8 7/108r1
5.6r1 C3H7 ? H2 7/102r1 1.4r1
H2
CH3
C3H7
C4H8
22
Elemfluxus számítása KINALC-al
c ATOMFLOW Fluxes of elements from
species to species are investigated c
The name(s) of elements are listed after the
keyword. c Usage ATOMFLOW ltelement1gt
ltelement2gt ... ATOMFLOW C H
ATOMFLOW
Fluxes of elements from species to
species Net fluxes of element H
absolute rel.
1 H2 gt H2O
6.843E-02 mole/(cm3 sec) 1.0000 2 H2
gt H 4.584E-02
mole/(cm3 sec) .6699 3 OH
gt H2O 4.034E-02 mole/(cm3
sec) .5895 4 H gt OH
3.360E-02 mole/(cm3 sec)
.4910 5 H gt HO2
2.370E-02 mole/(cm3 sec) .3463 6
OH gt H 2.302E-02
mole/(cm3 sec) .3364 7 HO2
gt OH 1.797E-02 mole/(cm3
sec) .2626 8 H2 gt OH
1.162E-02 mole/(cm3 sec)
.1699 9 H gt H2
6.346E-03 mole/(cm3 sec) .0927 10
H gt H2O 4.084E-03
mole/(cm3 sec) .0597 11 HO2
gt H2 3.334E-03 mole/(cm3
sec) .0487 12 HO2 gt H2O
2.689E-03 mole/(cm3 sec)
.0393 13 OH gt H2
1.049E-03 mole/(cm3 sec) .0153 14
H2O gt OH 7.891E-04
mole/(cm3 sec) .0115 15 OH
gt H2O2 7.617E-04 mole/(cm3
sec) .0111 16 H2O gt H2
7.010E-04 mole/(cm3 sec)
.0102 17 H2O2 gt H2O
5.367E-04 mole/(cm3 sec) .0078 18
H2O2 gt OH 2.131E-04
mole/(cm3 sec) .0031
23
KINALC ? FluxViewer

• KINALC eredményét jobb lenne ábrán (is) látni.
• ? FluxViewer JAVA program az elemfluxusok
  megjelenítésére

• Az anyagok cimkéi mozgathatók
• A nyilak vastagsága a
• log elemfluxusokkal arányos
• a nyilak vastagsága
• változtatható
• rajzok vagy mozifilm

24
C-fluxusok, ?1.0, metán-levego láng
T 815 K
25
C-fluxusok, ?1.0, metán-levego láng
T1155 K
26
C-fluxusok, ?1.0, metán-levego láng
T1500 K
27
C-fluxusok, ?1.0, metán-levego láng
T1805 K
28
C-fluxusok, ?1.0, metán-levego láng
T1865 K
29
C-fluxusok, ?1.0, metán-levego láng
T1915 K
30
Különbözo elem-fluxusok összehasonlítása metán-lev
ego robbanás (? 1.0, T1800 K)
C-fluxusok
31
Különbözo elem-fluxusok összehasonlítása metán-lev
ego robbanás (? 1.0, T1800 K)
O-fluxusok
32
Különbözo elem-fluxusok összehasonlítása metán-lev
ego robbanás (? 1.0, T1800 K)
H-fluxusok
33
3. Idoskála-analízis
34
3 Idoskála-analízis 3.1 Élettartamok és idoskálák
felezési ido Ennyi ido alatt csökken a felére
egy anyag koncentrációja, ha nem termelodik és az
összes többi anyag koncentrációja változatlan.
élettartam Ennyi ido alatt csökken az e-ed
részére egy anyag koncentrációja, ha nem
termelodik és az összes többi anyag
koncentrációja változatlan.
Egyetlen elsorendu reakció A ? P Élettartam
Felezési ido
Több elsorendu reakció (fotokémia,
fénygerjesztett részecske reagál) A ? P1 ?
P2 ? P3 Élettartam
35
Élettartam
légkörkémia gyökkoncentrációk kicsik, ezért
gyök gyök reakciólépések (pl. 2 CH3 ? C2H6)
hiányoznak a légkörkémiai mechanizmusokból ? Yi2
tagok nincsenek a diff. egyenletben Pi termelo
reakciólépések sebességének összege Li fogyasztó
reakciólépések sebességének összege Yi
koncentrációváltozási sebessége Yi élettartama
tetszoleges reakciómechanizmus Nincs korlátozás
a diff egyenlet polinom fokszámában Yi
élettartama ahol jii a Jacobi mátrix
foátlójának i-edik eleme.
36
Lassú változó - gyors változó
Ha egyetlen anyagfajta koncentrációját
változtatjuk meg ?yi-vel és a többi anyag
koncentrációja nem változik meg ennek hatására
rövid élettartamú anyagfajta a perturbáció
hatása gyorsan lecseng az eredeti és a
megzavart trajektória gyorsan közeledik
egymáshoz ? gyors változó
hosszú élettartamú anyagfajta a perturbáció
hatása lassan cseng le az eredeti és a
megzavart trajektória közel párhuzamosan halad ?
lassú változó
Következmények a gyors változók elfelejtik kezdo
értéküket a gyors változók értékét teljesen
megszabja a többi változó a gyors változók-lassú
változó besorolás független attól, hogy dYi
/d t mekkora
37
Lassú sokaságok dinamikus rendszerekben
A gyors módusok gyorsan relaxálódnak a
trajektória rázuhan egy altérre, majd ahhoz
nagyon közel fog haladni. 2D a trajektóriák egy
síkra zuhannak rá 1D a trajektóriák egy görbe
mentén helyezkednek el 0D elértük a stabil
egyensúlyi pontot
38
Lassú sokaságok
39
Lassú sokaság a koncentrációtérben
c2
1D lassú sokaság
c1
40
Jacobi-mátrix sajátérték-sajátvektor felbontása
sajátértékeket tartalmazó diagonális mátrix
(komplex sajátértékek!)
baloldali sajátvektorok (sorvektorok) mátrixa
jobboldali sajátvektorok (oszlopvektorok) mátrixa
Jelölje Wf azon sajátvektorok mátrixát,
amelyekhez tartozó Re(?) kicsi (negatív) szám
A bal- és jobboldali sajátvektorok ortonormáltak
Ennek következtében
41
Lassú sokaságok számítása
Legyen yM a sokaság felületén levo pont Ekkor
f(yM) a pont mozgási sebessége a sokaság
felületén. Wf rásimul a sokaság felületére,
ezekben az irányokban biztosan nem mozog a
pont f(yM) vektor ortogonális minden wf
vektorra
42
Módusok
Ha egyetlen anyagfajta koncentrációját
változtatjuk meg ?yi-vel és a többi anyag
koncentrációja nem változik meg ennek
hatására Csak akkor teljesül, ha a
megváltoztatott anyagfajta koncentrációja kicsi
és élettartama rövid.
Egyébként egy anyagfajta koncentrációjának
megváltoztatása a többi koncentrációváltozását is
kiváltja ? vizsgáljuk meg, mi van akkor, ha
egyszerre több anyagkoncentrációt változtatunk
meg ?Y-al
Diff. egyenlet (lineáris közelítéssel)
Diff. egyenlet megoldása, ha feltételezzük, hogy
a relaxációs ido alatt a Jacobi-mátrix
gyakorlatilag állandó maradt
?
ahol ? t - t1 a perturbáció óta eltelt ido
43
Módusok 2.
Bevezetünk egy új változókészletet z módusvektor
zj módus számítása
ahol wi a W mátrix i-edik sora
yi koncentráció visszanyerése
ahol vi a V mátrix i-edik oszlopa
A kinetikai diff. egyenlet transzformálva
Ha a ?y perturbáció w1 irányba történik, akkor
z1?0, de a többi j-re zj0.
FONTOS a Jacobi a koncentrációtérben helyrol
helyre változik, emiatt az y ? z transzformáció
is helyrol helyre változik!
44
Módusok 3.
Jelölje zf a gyors módusokat és yM a sokaságon
levo pontot.
Mivel
Csak az i-edik módusra
45
Módusok 4.
Ha csak a i-edik módus irányában perturbálunk,
akkor onnan egyexponenciális a visszatérés
Ennek t szerinti deriváltja
Az elozo oldal szerint
ahol ?zi a távolság a sokaságtól.
46
Dinamikus dimenzió számítása
wj irányban a rendszer állapotának távolsága a
lassú sokaságtól
Ha
? ez a módus (csaknem) a sokaságon van.
Tegyük fel, hogy nr darab ilyen sokaságra
relaxált módus van.
47
Stacionárius rendszer stabilitásvizsgálata
A dinamikus rendszert leíró diff. egyenlet
A stacionárius pontban
Stabilis stacionárius pont ha kimozdítjuk,
visszamegy

minden sajátértékének valós része negatív
Jacobi mátrix
Instabilis stacionárius pont ha kimozdítjuk,
eltér

Jacobi mátrix legalább egy sajátértékének valós
része pozitív
48
Mozgó rendszer stabilitásvizsgálata
A dinamikus rendszert leíró diff. egyenlet
Stabilis trajektória ha kimozdítjuk, visszamegy
az eredeti útvonalra
minden sajátértékének valós része negatív
Jacobi mátrix
Instabilis trajektória ha kimozdítjuk, eltér az
eredeti útvonaltól
Jacobi mátrix legalább egy sajátértékének valós
része pozitív
49
A hidrogén-levego adiabatikus robbanás dimenziója
pozitív Jacobi sajátértékek
n a változók száma a modellben nc a megorzött
tulajdonságok száma a Jacobi nulla
sajátértékeinek száma nr a relaxált módusok
száma Pillanatnyi dinamikus dimenzió
nDn?nc?nr.
50
4. Érzékenységanalízis
51
4 Érzékenység-analízis 4.1 lokális
érzékenység-analízis

• Az érzékenységanalízis matematikai módszerek egy
  családja.
• Azt vizsgálja, hogyan függ a modellek eredménye a
  paraméterek értékétol.

Lokális érzékenységanalízis paraméterek kis
megváltoztatásának hatása.
Lokális érzékenységi együttható számítása véges
differencia közelítéssel
t1 idopontban megváltoztatjuk a paramétert és t2
idopontban megnézzük a változtatás hatását.
52
Lokális érzékenységanalízis
Egy másik megközelítés Taylor-sor közelítés
Lokális érzékenységi együttható
Lokális érzékenységi mátrix
A lokális érzékenységek alapján becsülheto a
paraméterváltoztatás hatása
Egy paraméter megváltoztatása
Több paraméter megváltoztatása
53
Lokális érzékenységanalízis 2
Ezt az egyenletet deriváljuk kj szerint, az
eredmény
Ugyanez mátrix-vektor felírással
ahol
közvetett hatás
közvetlen hatás
54
Kezdeti koncentráció szerinti érzékenységek
Kezdeti koncentrációk szerinti érzékenységek
véges differenciával
Kinetikai differenciálegyenlet-rendszer
Ezt az egyenletet deriváljuk Yj(t1) szerint, az
eredmény
t1 idopontban megváltoztatjuk az j-edik változó
értékét és a t2 idopontban megnézzük a
változtatás hatását.
Green-függvény mátrix G
55
Lokális érzékenységek számítása
1. A véges differencia közelítés szerinti
számítás (brute force method a nyers eroszak
módszere)
?pj kicsi számábrázolás okozta hiba ?pj nagy
nemlinearitás okozta hiba
2. Direkt módszer (direct method) 2a. Csatolt
direkt módszer (coupled direct method) a
kinetikai és az érzékenységi egyenletek megoldása
együtt
Az együttes megoldás megismétlése minden
paraméterre
Mindenképpen felesleges számításokkal jár.
56
Lokális érzékenységek számítása 2.
2b. Szétcsatolt direkt módszer (decoupled direct
method) a kinetikai és az érzékenységi
egyenletek megoldása együtt
A fenti egyenletek Jacobi mátrixa azonos ? Jacobi
mátrix háromszög mátrix-szá alakítása ? Jacobi
mátrix alapján ?t lépéshossz választása ? Y
változóvektor új értéke az új idopontban ? a
korábban kiszámított háromszögmátrix
felhasználásával az érzékenységi
diffegyenlet megoldása ? újra lépés az idoben ?
? ? ?
gyors algoritmus - több száz paraméter esetén
is csak 2-3-szor több gépidobe kerül az összes
érzékenységi együttható kiszámítása -
érzékenységi együtthatók pontossága becsülheto
57
Lokális érzékenységek értelme
Egységnyi paraméterváltoztatás hány egységnyi
eredmény változást okoz? eredmény egysége /
paraméter egysége
Normált érzékenység
1 paraméterváltoztatás hány eredmény
változást okoz? dimenziómentes
Eddig egyetlen paraméter megváltoztatásának
hatását nézzük egyetlen modelleredményre
Érdekelhet minket egyetlen paraméter
megváltoztatásának hatása egyszerre
több modelleredményre
Bruttó érzékenység
58
PCAS az érzékenységi mátrix fokomponens analízise
Ha egyszerre több paramétert változtatunk és
több modeleredmény megváltozását figyeljük A
paraméterváltoztatás hatását egy célfüggvényen
mérjük le
59
PCAS az érzékenységi mátrix fokomponens
analízise 2.
A célfüggvény az alábbi alakkal közelítheto
ahol a tr idoponthoz tartozó normált
érzékenységi mátrix
60
PCAS az érzékenységi mátrix fokomponens
analízise 3.
Ez a kvadratikus alak egy (hiper)ellipszisthatáro
z meg 2D ellipszis3D rögbi-labda4D hiperrögbi-
labda
A (hiper)ellipszist egy másik egyenértéku
megadása- a tengelyek hossza- a tengelyek
iránya a paramétertérben
mátrix sajátérték-sajátvektor
felbontása ?i i-edik sajátérték i-edik
tengely hossza ha ?i kicsi abba az irányba
gyorsan emelkedik a célfüggvény i-edik
paramétercsoport hatásosui i-edik sajátértékhez
tartozó sajátvektor i-edik tengely iránya
61
PCAS az érzékenységi mátrix fokomponens
analízise 4.
A célfüggvény másik megadása
Ahol
A fokomponensnek nevezett transzformált
paramétervektor
Példa ?1 nagy u1 (0,707, 0,707) ?2 kicsi
u1 (-0,707, 0,707)
Megjegyzés a sajátvektorok egységvektorok 0,707
20,7072 1
62
PCAS az érzékenységi mátrix fokomponens
analízise 5.
Példa ?1 kicsi 1. tengely hosszú u1
(0,707, 0,707) ?2 nagy 2. tengely rövid u2
(-0,707, 0,707)
1. tengely hosszú ? u1 irányban megváltoztatva a
paramétereket sokáig nem változik a
célfüggvény ? ha ?2-?1 ln p2 - ln p1 ln
(p2/p1) állandó, akkor alig változik a
célfüggvény ? ha p2/p1 állandó, akkor alig
változik a célfüggvény
Tehát u (0,707, 0,707) sajátvektor azt jelenti,
hogy a megfelelo paraméterek hányadosát
állandónak tartva nem változik meg a szimuláció
eredménye
63
PCAS az érzékenységi mátrix fokomponens
analízise 6.
Példa 2 u1 ( 0,707, 0,707, 0) u2 (-0,707,
0,707, 0) u3 ( 0 , 0 , 1)
Egy ilyen PCAS eredményt adó rendszerben hogyan
kell paramétert becsülni?
Ha mindhárom paramétert engedjük változni, akkor
végtelen ideig szöszmötöl azon, hogy kicsit
növeli az 1. majd kicsit növeli a 2. paramétert.
Megoldás p1 1 legyen rögzített, p2 és p3
becsült paraméterek. p3 függetlenül becsülheto,
de a p2-re vonatkozó érték tulajdonképpen a p1/p2
hányados értékét adja meg, úgy is kell értelmezni!
64
Lokális érzékenységek felhasználása

• Modellek elemzése
• Paraméterperturbáció hatása a modellszámításra
• Kooperáló paraméterek azonosítása
• Modellek redukálása
• Hatástalan paraméterek azonosítása, és így sokkal
  kevesebb paramétert tartalmazó, ám elegendon
  pontos modellek eloállítása
• Lokális bizonytalanságanalízis
• A globális bizonytalanságanalízisnél kevésbé
  pontosan, de nagyon kis számításigénnyel ad
  eredményeket
• Paraméterbecslés
• A gradiensmódszerek mindig az érzékenységi
  együtthatók (rejtett) számításán alapulnak
• Hatásos paraméterek számának meghatározása
• Kísérlettervezés

65
4 Érzékenység-analízis 4.2 Globális
érzékenység-analízis
Lokális érzékenység analízis Információt ad
egy adott paraméterkészletnél Jól használható, ha
a vizsgált paramétertartományban nincs minoségi
változás
Globális érzékenység analízis Átvizsgálunk egy
(véges) paramétertartományt
A lokálishoz képest a globálishoz mindig sok
gépido kell. Megszerzett információ gépido
66
Globális érzékenység-analízis
A paraméterek bizonytalanságát valószínuségi
suruségfüggvénnyel (probability density
function, pdf) jellemezhetjük.
Globális érzékenységanalízis feladatai 1. A
paraméterek pdf-je következtében mi az eredmények
pdf-je? 2. Az eredmények szórásának mekkora
hányadát okozza egy-egy paraméter?

67
Morris-féle módszer
Rostáló módszerek elsodleges közelíto információ
gyorsan Nagy paraméterváltoztatások hatását
lehet vizsgálni

• minden paraméterhez megadjuk annak alsó és felso
  határát.
• n részre osztjuk minden paraméter intervallumát
• kisorsolunk egy véletlen paraméterkészletet
• minden további futás elott egyetlen paramétert
  változtatunk meg
• a futási eredményeket statisztikailag
  kiértékeljük
• nem használja fel a paraméterek val.
  suruségfüggvényét
• emiatt nem adja meg a megoldás suruségfüggvényeét
• közepes számítási igény

68
Morris-féle módszer 2
69
Monte Carlo módszer
Olyan, mint rulettezni Monte Carlo-ban -)
Kisorsolunk több ezer paraméterkészletet úgy,
hogy megfeleljen a paraméterek közös
pdf-jének. Elvégezzük a szimulációkat A futási
eredményeket feldolgozzuk, pl. a megoldás
hisztogramjának elkészítése a megoldás várható
értékének és szórásának számítása Sok
számítógépidot igényel. Nehéz kiszámítani
(elkülöníteni) az egyes paraméterek okozta hatást!
70
Monte Carlo módszer
100 db (?,?) pont ? és ? véletlen számok 0,1
egyenletes eloszlással
Csomósodás és üres foltok!
p1 hisztogramja
71
Monte Carlo módszer Latin hiperkocka
mintavétellel
100 db (?,?) pont latin hiperkocka mintavétellel
Sokkal egyenletesebb eloszlás!
p1 hisztogramja
72
Latin hiperkocka mintavétel










egyenletes eloszlás
normális eloszlás
Egyenlo valószínuségu sávokat (strata) jelölünk
ki. Egy-egy sávon belül véletlenszeruen jelölünk
ki pontot. Ha valahol van már pont, oda már nem
kerül (8 bástya probléma).
73
Fourier Amplitude Sensitivity Test (FAST) módszer
s változtatásával tekerünk nem összemérheto
frekvenciájú sin függvényeket
150 pont, ?s 0,1 x 0,5 (sin(17s),
sin(113s))0,5
74
Fourier Amplitude Sensitivity Test (FAST)
módszer 2.
Yi modell eredmény várható értéke E(Yi)
Ahol Yi kiszámítható hi kiértékelésével P pedig
p közös pdf-je
A pj paramétert az s skalár tekerésével
változtatjuk
Gj megfelelo megválasztásával tudjuk megkapni P-t.
?j a pj paraméterhez tartozó frekvencia Ezeknek
prímeknek kell lenniük. Ha -?ltslt? és N pontot
helyezünk el, akkor ezek a pontok a paramétertér
minden pontjához közel lesznek
75
Fourier Amplitude Sensitivity Test (FAST)
módszer 3.
Az N modelleredményt Fourier analízisnek vetjük
alá
Ahol ?2(Yi) az eredmény szórása, Ail és Bil a
Fourier együtthatók
Ha a Fourier együtthatókat az ?j frekvenciáknál
és felharmonikusaiknál számítjuk, akkor a j-edik
paraméter okozta parciális szórást kapjuk meg
76
Fourier Amplitude Sensitivity Test (FAST)
módszer 4.
A parciális szórás hányada
FAST lassú N 1,2 k2,5 50 paraméter 21000
szimuláció
A legfontosabb többletinformáció MC-vel
ellentétben felhasználjuk a szimulációk
sorrendjét is! (ez többlet gépidot nem okoz)
77
Érzékenységi indexek
Saltelli, A., Comput. Phys Commun., 145, 280(2002)
FAST egyik továbbfejlesztett változata. Yi
várható értéke
Yi szórásnégyzete
Yi szórásnégyzete, ha a pj paramétert rögzítjük
V(Yi?pj) Ennek várható értéke E(V(Yi?pj))
Yi szórásnégyzete, amit pj okoz
V(E(Yi?pj))V(Yi)- E(V(Yi?pj)) Az elsorendu
bizonytalansági index ( FAST parciális szórás)
78
Érzékenységi indexek 2.
Yi szórásnégyzete, amit pj és pk együttesen
okoz V(E(Yi?pj,pk)) Ebbol számítható a
másodrendu bizonytalansági index
Ez megadja pj és pk kölcsönhatását. Hasonlóan
számítható minden n-edrendu bizonytalansági
index.
Legyen három paraméterünk a, b, c a teljes
indexe
Ha a j index semmi mással nem korrelál
j index kölcsönhatásai
79
Érzékenységi indexek 3.

• globális módszer
• ez is álvéletlen-számokat használ, hogy az
  integrálok könnyen kiszámolhatók legyenek
• paraméterek elsodleges és magasabb rendu hatása
• a totális hatást is számolja
• a paraméterek pdf-jét is tekintetbe veszi
• nagyon számításigényes
• (50 paraméterre kb. 25000 futás)

80
A módszerek összehasonlítása
local Morris MC LHS érz. index
input szórása ? ? ? ?
input pdf ? ? ? ?
output pdf ? ? ? ?
output szórása ?(lineáris) ? ? ?(torzított)
Gépido igény? ? 1 2110 3000 16280
egyedi hozzájárulások ? (lineáris) ?(csak kvalitatív) ? ?
globális? ? ? ? ?
Info nemlinearitásról ? ? (csak kvalitatív) ? ?
81
5. Érzékenységi függvények hasonlósága
82
Mi az, hogy függvények hasonlósága?
Az érzékenységi függvények között nincs a priori
kapcsolat, általában így néznek ki.
Többen is észrevették, hogy bizonyos
rendszerekben az érzékenységi függvények nagyfokú
szabályosságot mutatnak.
83
1. szabályosságAz érzékenységi együtthatók
lokális hasonlósága
pirosat a pirossal, zöldet a zölddel osztjuk
84
1. szabályosságAz érzékenységi együtthatók
lokális hasonlósága
85
2. szabályosságSkálaviszony törvény
86
2. szabályosság a skálaviszony-törvény
87
3. szabályosságAz érzékenységi együtthatók
globális hasonlósága
a pirosat a zölddel osztjuk
88
3. szabályosságAz érzékenységi együtthatók
globális hasonlósága
állandó a független változó széles tartományában.
89
A lokális hasonlóság és a skálaviszony eredete1D
sokaságon geometriai bizonyítás

• Szükséges feltételek
• a rendszer 1D sokaságon mozogjon
• az 1D sokaság kevéssé mozduljon el a paraméterek
  megváltoztatásakor

90
Mozgás nem 1D sokaságon
Nem feltétlenül teljesül, hogy a
paraméterperturbáció hatása párhuzamos a
trajektóriával ? nem feltétlenül lesz lokális
hasonlóság, illetve skálaviszony
Továbbá sokaság dimenziója ? az érzékenységi
mátrix rangja
91
A lokális hasonlóság és a skálaviszony eredete 1D
sokaságon bizonyítás differenciálszámítással
Az 1D sokaság egyenlete Y1 a paraméterezo
változó Fi az i-edik változó értéke z a
független változó (pl. ido)
Deriválás z szerint
Deriválás pj szerint
A skálaviszony egyenlete
92
A sokaság dimenziója ? az S mátrix rangja
bizonyítás differenciálszámítással
Az nD sokaság egyenlete Y1, Y2, ..., Yn a
paraméterezo változók Fi az i-edik változó
értéke z a független változó (pl. ido)
Deriválás pj szerint
Az S mátrix rangja legfeljebb n
93
A globális hasonlóság eredete bizonyítás 1
a rendszert leíró k.é.p.
az érzékenységi k.é.p.
Érzékenység számítás GFM-en keresztül
A második tag a pszeudo-homogén tartományban
elhanyagolható
94
A globális hasonlóság eredete bizonyítás 2
felhasználjuk a lokális hasonlóságot
a jobb oldal paraméter-független!
ez a globális hasonlóság
95
Igaz-e, hogy a homérséklet az egyetlen
domináns változó égéseknél?
Ha az egyik változó domináns
96
H2levego elegyek robbanása
Adiabatikus robbanás ? rögzített
homérsékletprofil (követi az adiabatikus
profilt)
97
H2levego elegyek lamináris lángja
Adiabatikus láng ? rögzített
homérsékletprofil (követi az adiabatikus
profilt)
98
6 különbözo égési rendszer együttes
koncentrációhomérséklet görbéi
homérséklet robbanások és lángok közös nevezoje
Robbanás Szabadon terjedo láng
Égofejhez rögzített láng --------
99
Érzékenységi együtthatók a homérséklet
függvényében égések esetén általában ilyenek
100
DE az adiabatikus hidrogénlevego elegy
robbanása esetén ...
nagyon szépek!
101
DE az adiabatikus hidrogénlevego elegy
robbanása esetén ...
pedig ...
H2 O OH H H2 OH H2O H H2O H H2
OH O2 H M HO2 M O2 H H2O HO2
H2O O2 H OH O OH O O2 H H2O2 H
HO2 H2 H2O2 H OH H2O H2O2 O OH
HO2 H2O2 OH H2O HO2 H2O2 (M) OH OH
(M) OH OH ( M) H2O2 ( M) H H M H2
M H H H2 H2 H2 H O M OH M H
OH M H2O M H HO2 H2 O2 H HO2 OH
OH H HO2 H2O O O HO2 O2 OH OH
OH O H2O OH HO2 H2O O2 HO2 HO2
H2O2 O2 HO2 HO2 H2O2 O2
9 anyagfajta 46 reakciója
102
A hidrogén adiabatikus robbanás dimenziója
pozitív Jacobi sajátértékek
103
Adiabatikus robbanás
Lokális hasonlóság VAN
Skálaviszony törvény VAN
Globális hasonlóság VAN
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
104
Rögzített homérsékletprofilú robbanás
Lokális hasonlóság VAN egyes paraméterekre
Skálaviszony törvény NINCS
Globális hasonlóság NINCS
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
105
Adiabatikus égofej stabilizált láng
Lokális hasonlóság VAN egyes paraméterekre
Skálaviszony törvény VAN egyes paraméterekre
Globális hasonlóság VAN egyes paraméterekre
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
106
Rögzített homérséklet profilú égofej stabilizált
láng
Lokális hasonlóság NINCS
Skálaviszony törvény NINCS
Globális hasonlóság NINCS
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
107
Adiabatikus szabadon terjedo láng
Lokális hasonlóság NINCS
Skálaviszony törvény NINCS
Globális hasonlóság NINCS
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
108
Rögzített homérséklet profilú szabadon terjedo
láng
Lokális hasonlóság NINCS
Skálaviszony törvény NINCS
Globális hasonlóság NINCS
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
109
Rögzített homérsékletprofil
Ha homérsékletprofil rögzített, akkor is elromlik
a hasonlóság, mert másik trajektóriára
kényszerítjük a rendszert.
110
A lokális hasonlóság, mint az érzékenységek
korrelációja
111
Pókháló ábrák

• Korrelációt mutatja
• Egyszerre több vektor hasonlítható össze

a b c
d
max.
0
min.
112
Pókháló ábrák

• Korrelációt mutatja
• Egyszerre több vektor hasonlítható össze

a b c d
113
A pókháló-ábrák elonye és hátránya

• Elonye
• Minden érzékenységi vektor egyenrangú, nem kell
  kijelölni vonatkoztatási vektort
• Látványosan megjelennek a korrelációk

• Hátránya
• Minden pókháló-ábra egy adott idoponthoz
  tartozik, az idobeni fejlodést nem mutatja meg
• Nem számszerusíti a korreláció mértékét

A korreláció-függvények ábrázolásának nincsenek
meg ezek a hátrányai. (De nincsenek meg a fenti
elonyei sem.)
114
H2levego adiabatikus robbanás pókháló ábrái
? 1.0
T 1935 K
115
H2levego adiabatikus robbanás pókháló ábrái
? 1.0
T 1940 K
116
H2levego adiabatikus robbanás pókháló ábrái
? 1.0
T 1945 K
117
A korrelációfüggvény
A matematikai statisztikában megszokott
korreláció-függvény Tartalmaz kompenzációt az
elemek átlagának különbségére és nyújtására

• Két érzékenységi vektor lokális hasonlóságának
  vizsgálatakor csak
• a nyújtás okozta eltérést kell vizsgálni
• a fenti képleteket egyszerusíthetjük úgy, hogy
  kivesszük azeltolást az átlaggal ? a maradék
  tulajdonképpen két egységnyi hosszranormált
  vektor skalárszorzata

A következo ábrákon a vonatkoztatási vektor
mindig a vízhez tartozó érzékenységi vektor.
118
Korreláció-függvények
? 1.0
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
119
Korreláció-függvények
? 1.0
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
120
Korreláció-függvények
? 1.0
A robbanás R robbanás A égofej láng R égofej
láng A szabad láng R szabad láng
121
Miért fontos a globális hasonlóság?
k1
k2
X
megoldás
122
3. szabályosságAz érzékenységi együtthatók
globális hasonlósága
a pirosat a zölddel osztjuk
123
Miért fontos a globális hasonlóság?

• ??? metán robbanás (eredeti kémiai mechanizmus)
• - - - - 4 paramétert elhangultunk 50-kal
• . egy 5. paramétert hangoltunk 9.875-kal

CO2 koncentráció
124
Miért fontos a globális hasonlóság?

• ??? metán robbanás (eredeti kémiai mechanizmus)
• - - - - 4 paramétert elhangultunk 50-kal
• . egy 5. paramétert hangoltunk 9.875-kal

OH koncentráció
CO2 koncentráció
CH3 koncentráció
CH koncentráció
AZ EREDETI GÖRBÉKET VISSZAKAPTUK MINDEN
VÁLTOZÓRA és MINDEN IDOPONTBAN
125
A globális hasonlóság jelentosége
Empirikus modellek Illeszthetok egyetlen
tetszoleges (hatékony) paraméterrel (csak
egyetlen egyet hangoltam a 35 felhasznált
paraméterbol) Fizikai modellek A használt
paraméterek jósága ellenorizheto összetett
modellekkel, DE a fizikai modellhez illesztett
paraméternek általában nincs fizikai értelme.
(lásd k meghatározása komplex kinetikai
rendszerekbol!!!) Sejtek önszabályozása minden
pontban minden anyagra külön szabályozó
rendszer VAGY globális hasonlóság esetén egyetlen
paraméter változtatásával minden pontban
egyszerre visszaállítható az összes változó
optimális értéke Egy új eszköz gyógyszerek
tervezésére Javítási módszerek jelenleg ami
elromlik, azt megjavítjuk. Biológiai rendszer
modellje ? megtalálni azokat a paramétereket,
amelyek globálisan hasonlók. A globálisan
hasonló paraméterek csoportján belül tetszoleges
paraméterrel meg tudjuk javítani a muködést.
126
Vajon a globális hasonlóság a dinamikai
rendszerek egy általános tulajdonsága?
A legtöbb fizikai és kémiai modell olyan
folyamatokat ír le, amelyek nagyon eltéro
idoskálán mozognak ? vannak lassú sokaságok ?
lokális hasonlóság
A dinamikus rendszerek az egyensúlyi ponthoz
közel 1D sokaságon vannak ? 1D lassú sokaság ?
skálaviszony
Lokális hasonlóság és az ODE pszeudo-homogenitás
a ? globális hasonlóság pszeudo-homogenitás ?
autokatalitikus folyamatok
127
1. rész vége