1
Model vznícení
zpracovala Hana Tužilová podle Andrew Fowlera
0.1 Model vznícení
Asi každý z nás už někdy zapaloval zápalku. Proč se ale zapálka při škrtnutí o krabičku
vznítí? Odpověď na tuto otázku budeme hledat pomocí vhodného modelu.
K tomu, aby se zápalka vznítila, je potřeba dodat určitou aktivační energii E. Poté se
rozběhne exotermická chemická reakce (reakce, během které se uvolňuje teplo). Množství
uvolněného tepla je úměrné rychlosti reakce a rychlost samotné reakce vzrůstá s rostoucí
teplotou (když zapálíme dřevo v ohništi, zpočátku oheň hoří pomalu a teprve postupně se
jeho intenzita zvětšuje). Teplo uvolňované během hoření je popsáno Arrheniovým vztahem
A exp −
E
RT
,
kde E > 0 je aktivační energie, R > 0 univerzální plynová konstanta, T teplota zápalky
v kelvinech a A > 0 konstanta (frekvenční faktor, ve skutečnosti závislý na koncentraci
reaktantů).
Současně dochází k tepelné výměně mezi zápalkou a okolím. Je-li teplota okolí T0, potom
podle Newtonova zákona ochlazování je tato výměna popsána vztahem
−k(T − T0)
s koeficientem ochlazování k > 0.
Jednoduchý model vývoje teploty zápalky potom vypadá následovně
c ·
dT
dt
= −k(T − T0) + A exp −
E
RT
, (1)
kde c je odpovídající měrná tepelná kapacita (množství tepla potřebné k ohřátí 1 kg látky
o 1 K).
Diferenciální rovnici (1) neumíme řešit analyticky, její řešení bychom mohli najít přibližně
prostředky numerických metod. Pro nás však konkrétní řešení této rovnice není stěžejní.
Na problém se budeme dívat graficky a model prozkoumáme prostřednictvím nástrojů
kvalitativní teorie diferenciálních rovnic.
Rovnovážné body jsou řešení rovnice
k(T − T0) = A exp −
E
RT
. (2)
Funkce na pravé, resp. levé, straně mají obvykle průběh znázorněný na obrázku1
1. Vidíme,
1
Protože u nás teplotu měříme ve ◦
C, pro rychlejší zorientování čtenáře je teplota ve všech zobrazených
grafech uvedena ve ◦
C.
2
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
−0.2
0
0.2
0.4
0.6
0.8
1
T
0
Obrázek 1: Funkce k(T − T0)
a A exp {−E/R(T + Tm)} pro hodnoty
A = 1, E = 20000, R = 8.3, k =
10−4
, T0 = 15◦
C, Tm = 273.
0 100 200 300 400 500 600 700 800 900 1000
0
1000
2000
3000
4000
5000
6000
T
0
∆
T
−
T
+
Obrázek 2: Křivka různých rovnovážných
bodů pro ∆ v závislosti na T0. Parametry
jsou stejné jako u obrázku 1, jen E =
35000.
že rovnice má ve zobrazeném případě tři rovnovážné body: dva stabilní (horní a spodní)
a jeden nestabilní. Dolní rovnováha odpovídá klidovému stavu (zápalka leží v krabičce)
a horní rovnováha tomu, když zápalka hoří. pro jinou volbu parametrů může samozřejmě
nastat situace, kdy bude existovat jen horní (např. pro velké hodnoty T0) nebo jen dolní
rovnováha (pro malé hodnoty T0).
Označíme-li ∆ := T − T0, rovnici (2) můžeme psát jako
f(∆, T0) := −k∆ + A exp −
E
R(∆ + T0)
= 0. (3)
Předmětem našeho zájmu bude vývoj rozdílu teploty zápalky a okolí ∆ při změně okolní
teploty T0. Řešení rovnice (3) vzhledem k neznámé ∆ neumíme vyjádřit explicitně, je však
zajímavé pozorovat, jak se bude řešení chovat, budeme-li jej hledat numericky. Pro numerický
výpočet řešení (3) si zvolíme vhodné hodnoty parametrů A, E, R, k (např. A =
1, E = 35000, R = 8.3, k = 0.001), a vybereme si nějaký výpočetní program – zde uvedeme,
jak by se takový výpočet provedl v programu MATLAB.
Nejdříve si vytvoříme m-soubor s funkcí f
function f=zadani(Delta,T0,E,A,R,k)
f=-k*Delta+A*exp(-E/R/(Delta+T0+273));
a budeme pro zástupce hodnot E, A, R, k hledat její řešení pro různá T0. Pro numerické
řešení je třeba specifikovat nějaký odhad hodnoty, kterou hledáme. My pro první výpočet
počáteční odhad ”střelíme od pasu” a poté již budeme za počáteční odhad brát vždy předchozí
stav. Kvůli názornosti si necháme vykreslit řešení jak pro rostoucí, tak pro klesající
teplotu T0.
3
>> E=35000; A=1; R=8.3; k=0.0001;
>> sol1(1)=fzero(@(Delta) zadani(Delta,1,E,A,R,k),4000);
>> for T0=2:1000
sol1(T0)=fzero(@(Delta) zadani(Delta,T0,E,A,R,k),sol1(T0-1)); end;
>> x=1:1000;
>> plot(x,sol1,’.r’)
>> sol2(1000)=fzero(@(Delta) zadani(Delta,1000,E,A,R,k),4000);
>> for T0=999:-1:1
sol2(T0)=fzero(@(Delta) zadani(Delta,T0,E,A,R,k),sol2(T0+1)); end;
>> plot(x,sol2,’.r’)
Výstup, zobrazený na obrázcích 3 a 4, můžeme interpretovat následovně. Řekněme, že
0 100 200 300 400 500 600 700 800 900 1000
0
1000
2000
3000
4000
5000
6000
T
0
∆
Obrázek 3: Řešení rovnice f(∆, T0) = 0
pro rostoucí teplotu T0.
0 100 200 300 400 500 600 700 800 900 1000
0
1000
2000
3000
4000
5000
6000
T
0
∆
Obrázek 4: Řešení rovnice f(∆, T0) = 0
pro klesající teplotu T0.
bychom dali velmi dlouhou zápalku do studené trouby a zvolna zvyšovali její teplotu T0.
Teplota zápalky pak nejprve pomalu narůstá, až se při teplotě, kterou označíme jako T+,
zápalka vznítí a dojde k teplotnímu skoku. Naopak budeme-li teplotu T0 postupně od vysokých
hodnot snižovat, bude teplota zápalky pomalu klesat, až při určité teplotě, označme
ji T−, zápalka zhasne a dojde opět k teplotnímu skoku.
Křivka všech rovnovážných bodů (i těch nestabilních) – řešení rovnice (3) – ovšem vypadá
jako na obrázku 2. Snadno se o tom můžeme přesvědčit. Z rovnice (3) vyjádříme T0 jako
funkci ∆
T0(∆) = −
E
R ln k∆
A
− ∆ (4)
a při vykreslování výsledku přehodíme vstupní a výstupní hodnoty, tzn. nebudeme vykreslovat
T0 jako funkci ∆, ale ∆ v závislosti na T0. V MATLABu pak stačí zadat sekvenci
příkazů
4
>> y=1:0.01:5000;
>> x=-E/R./log(k*y/A)-y;
>> plot(x-273,y,’r.’)
Nabízí se jistě možnost podrobně vyšetřit průběh funkce T0(∆) pomocí prvních a druhých
derivací, to nicméně vede k dalším rovnicím, které nelze řešit jinak než numericky.
Určeme ještě, při jakých hodnotách vlastně dochází ke vznícení a ke zhasnutí zápalky.
Jedná se o body, ve kterých je derivace (4) rovna nule, podrobněji hledáme ∆ splňující
E
R∆
·
1
ln2 k∆
A
− 1 = 0. (5)
Opět pomocí MATLABu bychom pro zde zvolené parametry dostali výsledek ∆+ ≈ 418◦
C
a ∆− ≈ 3136◦
C, odkud již snadno spočítáme, že T+ ≈ 637◦
C a T− ≈ 227◦
C.
Z obrázku 2 je patrné, že se v tomto modelu projevuje jev hystereze ve smyslu nepřevratitelnosti;
cesta není oblouk, ale uzavřená křivka.
Vraťme se k původní otázce – proč zápalka začne hořet? Když škrtneme zápalkou, hlavička
zápalky se třením zahřeje z pokojové teploty na teplotu takovou, že ∆ = T −T0 bude větší
nebo rovno ∆+ (hodnota ∆ příslušná T+), zápalka se vznítí a ∆ vzroste až na hodnotu
v horním stabilním rovnovážném bodě. Třením tedy dodáváme aktivační energii potřebnou
pro rozběhnutí reakce.
Obrázek 2 rovněž vysvětluje, proč je tak obtížné zapálit mokrou zápalku (mokrá zápalka
má vyšší měrnou tepelnou kapacitu, což se projeví menší změnou teploty zápalky T při
dodání stejného množství energie; pro zapálení mokré zápalky je tedy třeba dodat větší
energii) a proč se zápalka sama vznítí, pokud ji dáme dostatečně blízko k hořící svíčce.