DIFFERENCIAEGYENLETEK, MINT A MODELLEZÉS ESZKÖZEI AZ ISKOLAI MATEMATIKÁBAN KOVÁCS ZOLTÁN
1. Bevezetés A természeti jelenségeket sokszor differenciálegyenletekkel lehet leírni: a vizsgált mennyiség például az idő függvénye, s erre a függvényre valamilyen differenciálegyenletet lehet felírni. Ezek közelítő megoldása differenciaegyenletek segítségével a középiskolás matematikai ismereteket nem haladja meg, így lehetővé válik a való életből vett bonyolultabb problémák matematikai modellezése. A differenciaegyenletek megoldásában, pontosabban a vizsgált mennyiség diszkrét időpontokban vett sorozatának előállítására, használhatunk akár komputeralgebrai eszközöket, akár táblázatkezelő alkalmazásokat. Mindkét esetben lehetőségünk van a grafikus reprezentációra is. Ebben a munkában a Maxima komputeralgebrai rendszert használom, wxMaxima 0.8.2 interfésszel. A teljes munkalapok elérhetők a honlapomon: zeus.nyf.hu/~kovacsz/maxima. A differenciaegyenletek középiskolás alkalmazása már a hetvenes években felvetődött az irodalomban. Lothar Berg [2] munkája 1982-ben magyarul is megjelent (Másodrendű differenciaegyenletek, Középiskolai szakköri füzet, Tankönyvkiadó, 1982), de a szerzőnek több német nyelvű munkája is megelőzte ezt (pl. [1]). Ezekben a munkákban főleg fizikai alkalmazások szerepeltek. A technológia alkalmazása új lendületet adhat ennek a lehetőségnek ([3]). A matematik tanítása során sokszor felteszik a tanárnak a kérdést: „ezt minek tanuljuk”? A való életből vett problémák tárgyalása a matematika tanításának egyik kulcskérdése: az igazi motivációt sok tanuló számára ez jelenti, tartalmat ad annak az idézetnek, amelyet sok szaktanterem falán láttam kiírva: „A természet a matematika nyelvén szól hozzánk”. 2. Kezdetiérték probléma numerikus megoldása A következő kezdetiérték problémát akarjuk numerikusan megoldani: (1)
y′ = f (t, y),
y(t0 ) = y0
a [t0 , b] intervallumon. Numerikus megoldáson azt értjük, hogy a megoldást véges sok pontban approximáljuk. Ezek a pontok pl. (2)
(t0 , y0 ), (t0 + ∆t, y1 ), (t0 + 2∆t, y2 ), . . . (t0 + n∆t, yn ),
ahol
b − t0 . n Az y′(tk ) derivált értéket a következő differenciával helyettesítjük: yk+1 − yk y′(tk ) ≈ ∆t ∆t =
1
2
KOVÁCS ZOLTÁN
Ami azt jelenti, hogy (1) helyett áttérünk a (3)
yk+1 = yk + ∆tf,
y(t0 ) = y0
differenciaegyenletre. A kritikus pont f argumentuma (3) jobb oldalán. Használhatjuk az f (tk , yk ) (első pont módszer ), f (tk+1 , yk+1 ) (hátsó pont módszer ), értékeket, ezek átlagát (szimmetrikus módszer ) vagy akár a [tk , tk+1 ] intervallumon vett több értéknek az átlagát is. Az első pont módszer szerint (4)
yk+1 = yk + ∆tf (tk , yk ),
tk = t0 + k∆t,
a hátsó pont módszer szerint (5)
yk+1 = yk + ∆tf (tk+1 , yk+1 ),
tk+1 = t0 + (k + 1)∆t,
míg a szimmetrikus módszer szerint (6)
1 yk+1 = yk + ∆t (f (tk , yk ) + f (tk+1 , yk+1 )) . 2
A továbbiakban a hátsó pont módszerrel nem foglalkozunk. Az első pont módszer egy egyszerű rekurziót határoz meg, azaz a kezdeti értékek alapján a keresett (2) sorozatot elő tudjuk állítani. A szimmetrikus módszer implicit egyenletet ad yk+1 re. Ha ennek megoldása elemi függvényekkel nehéz vagy nem lehetséges, akkor yk+1 -t az első pont módszer szerint közelítjük: (7)
1 yk+1 = yk + ∆t (f (tk+1 , yk + ∆tf (tk , yk )) + f (tk , yk )) . 2
A differenciálegyenletek numerikus megoldásának fentebb vázolt módszerét Eulermódszernek is nevezik. A következőekben példákkal illusztrálom a módszert. 3. Korlátozott növekedés Jelölje y(t) egy populáció egyedeinek számát. Semmilyen populáció nem növekedhet korlátlanul, korlátozást jelenthet az élelem hiánya, a szűk territórium, stb. Tegyük fel, hogy egy territórium B egyedszámot tud fenntartani és a populáció növekedése a B − y különbséggel arányos. (Minél inkább megközelíti a populáció a maximális egyedszámot, annál kisebb a növekedés.) Ezt a növekedési folyamatot a ∆y = K(B − y), ∆t
y(0) = y0
differencia egyenlet írja le, ahol K > 0 statisztikusan meghatározott arányossági tényező. Az egyenletet az első pont módszerrel megoldva: yk+1 = yk + ∆t(B − yk ) A sorozat előállításásnak Maxima kódja K = 0, 4, B = 2500 paraméterértékekkel, x0 = 150 induló egyedszámmal, ∆t = 0, 1 választása mellett 200 diszkrét időpontban: K:0.40$ B:2500$ jobb_oldal(y):=float(K*(B-y))$ dt:0.1$ N:200$
DIFFERENCIAEGYENLETEK. . .
3
x0=150 x0=1000
2500
2000
x
1500
1000
500
0 0
5
10
15
20
25
t
1. ábra. Korlátozott növekedés különböző induló egyedszámnál
y0:150$ list:[[0,y0]]$ y:y0$ for i:1 thru N-1 do(y:y+dt*jobb_oldal(y), list:append(list, [[(i+1)*dt,y]]))$ A grafikon kirajzolása: plot2d([’discrete,list], [style,[points,3,1,3]], [xlabel,"t"],[ylabel,"y"]); Az 1. ábrán a grafikont különböző induló egyedszámmal készítettem el (150 és 1500). Talán meglepő, hogy tízszeres induló populációszám alig jelent különbséget a territórium teljes benépesítésének idejében. 4. Lineáris harmonikus oszcillátor Az egydimenziós mozgásegyenlet az egység tömegű anyagi pontra: (8)
x′ = v,
v′ = F
ahol v a sebesség, F az egység tömegű anyagi pontra ható erők összege (m = 1). Ha az anyagi pont egyetlen rugó hatása alatt mozog, akkor a rugó erőtörvényét az F = −K ·x Hook-törvény írja le. Az egyenlet megoldására a szimmetrikus módszert
4
KOVÁCS ZOLTÁN
alkalmazzuk 1 xk+1 = xk + ∆t (vk + vk+1 ) 2 1 vk+1 = vk + ∆t (Fk + Fk+1 ) . 2
(9)
Az implicit egyenletet elkerülendő vk+1 -t az első pont módszerrel közelítjük:
(10)
1 1 xk+1 = xk + ∆t vk + ∆t (vk + ∆tFk ) 2 2 1 = xk + ∆tvk + ∆t2 Fk . 2
vk kiküszöböléséhez alkalmazzunk index transzformációt: (11)
1 xk+2 = xk+1 + ∆tvk+1 + ∆t2 Fk+1 . 2
Kivonva egymásból az (10) és (11) egyenleteket: 1 xk+2 = 2xk+1 − xk + ∆t (vk+1 − vk ) + ∆t2 (Fk+1 − Fk ) . 2 Helyettesítsük be vk+1 − vk az (9) egyenletből: xk+2 = 2xk+1 − xk + ∆t2 Fk+1 . A Hook-törvényt beírva: (12)
xk+2 = −xk + (2 − ∆t2 K)xk+1 .
Az hxk i sorozatra ismét egy rekurzív szabályt kaptunk, de most a sorozat egy elemét a két előtte lévő elmből kapjuk meg. Így x0 mellett szükségünk van x1 értékére is, melyet a (10) egyenletből kapjuk 1 x1 = x0 + ∆tv0 − ∆t2 Kx0 . 2 Részlet a Maxima kódból: force(x):=float(-K*x)$ nx1:float(x0+v0*dt+1/2*dt^2*force(x0))$ listb:[[0,x0],[dt,nx1]]$ for i:1 thru N do( listb:append(listb, [[(i+1)*dt,2*listb[i+1][2]-listb[i][2]+ force(listb[i+1][2])*dt^2]]))$ A 4. ábrán ugyanazon kezdeti feltételek mellett az itt nem részletezett (de az olvasó által könnyen levezethető) első pont módszer és a szimmetrikus módszer eredménye látható. Az első pont módszer mint modell nyilvánvalóan nem elfogadható, mert az enegiamegmaradás törvénye nem teljesül. Megjegyzendő, hogy a Hook-törvény esetén (9)-ből az hxk i sorozat az első pont módszer alkalmazása nélkül is kifejezhető, mert egy egyszerű lineáris egyenletrendszert kapunk.
DIFFERENCIAEGYENLETEK. . .
6
5
elso pont szimmetrikus
4
2
x
0
-2
-4
-6
-8 0
2
4
6 t
8
10
12
2. ábra. Harmonikus lineáris oszcillátor.
Hivatkozások [1] Lothar Berg. Mechanik ohne Differentialgleichungen – verwirklicht durch diskrete Modelle. Wissenschaft und Fortschritt, 28:146–149, 1978. [2] Lothar Berg. Differenzengleichungen zweiter Ordnung mit Anwendungen. VEB Deutscher Verlag der Wissenschaften, Berlin, 1979. [3] Zoltán Kovács. Difference equations as a modelling tool in school mathematics. In Proceedings of the IMEM 2009 Congress, Spisska Kapitula, ISBN 978-80-8084-471-4, pages 586–594, Catholic University in Ruzemberok, 2009. online: http://zeus.nyf.hu/˜kovacsz/shortpub/index.html. Nyíregyházi Főiskola, Matematika és Informatika Intézet, 4401 Nyíregyáza, Pf. 166 E-mail address:
[email protected]
