Lijnmethode

Uit Wikipedia, de vrije encyclopedie
Ga naar: navigatie, zoeken

De Lijnmethode is een numerieke methode om een parabolische partiële differentiaalvergelijking, of een aantal gekoppelde parabolische partiële differentiaalvergelijkingen op te lossen. De ruimtelijke coördinaten worden gediscretiseerd door middel van eindige differenties waardoor een stelsel van gekoppelde gewone differentiaalvergelijkingen ontstaat dat met de bekende numerieke methoden zoals de Runge-Kuttamethoden, kan worden opgelost. Als beginvoorwaarden moet de oplossing op t = 0 gegeven worden, en als randvoorwaarden de waarden op elk van de ruimtelijke randen.

Algemeen principe[bewerken]

De beschrijving die in deze paragraaf wordt gegeven is specifiek voor twee ruimtelijke dimensies. De gevallen van één of drie ruimtelijke dimensies kunnen op eenvoudige manier worden gevonden, respectievelijk door de y-afhankelijkheid te schrappen, of door een derde z-afhankelijkheid toe te voegen. Een parabolische partiële differentiaalvergelijking in twee ruimtelijke dimensies heeft als standaardvorm:

 \frac{\partial u}{\partial t} = k \, ( \frac{\partial^2 u}{\partial x^2} \, + \, \frac{\partial^2 u}{\partial y^2} ) \!

Deze vergelijking staat ook bekend als de diffusievergelijking of de vergelijking voor warmtetransport. In de Lijnmethode (in het engels Method of Lines) worden de ruimtelijke assen gediscretiseerd tot een verzameling netpunten, gekenmerkt door een stapgrootte h die hier in beide dimensies gelijk verondersteld wordt, zoals meestal het geval is in praktische toepassingen:

(x_i,y_j) \, \, \, \, \!

met

 \, \, \, x_i=i.h \quad (i =0..N) \!

en

 \, \, \, y_j=j.h \quad (j =0..M) \!

Door gebruik te maken van deze netpunten kunnen de afgeleiden gediscretiseerd worden, gebruik makend van eindige differenties. Een tweede orde afgeleide wordt daarbij omgezet in de centrale eindige differentie

 \frac{\partial^2 u}{\partial x^2}(t,x_i) \, \approx \, \frac{u(t,x_{i-1})-2u(t,x_i)+u(t,x_{i+1})}{h^2}

Toegepast op de diffusievergelijking met twee ruimtelijke dimensies wordt dit:

 \frac{\partial u}{\partial t} \, \approx \, k \, \frac{u(t,x_{i-1},y_j)+u(t,x_i,y_{j-1})-4u(t,x_i,y_j)+u(t,x_i,y_{j+1})+u(t,x_{i+1},y_j)}{h^2}

of korter:

 \frac{\partial u_{i,j}}{\partial t} \, = \, \frac{k}{h^2} \, (u_{i-1,j} \, + u_{i,j-1} \, -  4u_{i,j} \, +  u_{i+1,j} \, +  u_{i,j+1} )

waarbij u_{i,j} het verloop in de tijd geeft van de onbekende functie u, in het netpunt bepaald door de indices i en j. Op die manier ontstaat in elk netpunt (i,j) een gewone differentiaalvergelijking, die gekoppeld is aan vier naburige differentiaalvergelijkingen (twee per ruimtelijke dimensie). Het stelsel bevat, mits de nodige begin- en randvoorwaarden gegeven worden (zie verder), in totaal (N-1)x(M-1) gekoppelde gewone differentiaalvergelijkingen. Deze kunnen dan met de voor dit soort vergelijkingen beschikbare methoden numeriek worden opgelost, zoals bijvoorbeeld de bekende Runge-Kuttamethode. Eén enkele partiële differentiaalvergelijking wordt dus omgezet in een groot aantal gewone eerste orde differentiaalvergelijkingen.

Begin- en randvoorwaarden[bewerken]

Als beginvoorwaarde, dient de waarde van u(t,x,y) gegeven te worden in elk netpunt. De benodigde beginvoorwaarde voor de gewone differentiaalvergelijking in netpunt (i,j) is dus :

u_{ij}(0) =  u(t=0,x=i.h,y=j.h) \!

Daarnaast moeten de nodige randvoorwaarden gegeven worden. Deze kunnen verschillende vormen aannemen.

  • Bij de zogenaamde Dirichlet-randvoorwaarden wordt gewoon de gewenste waarde van de onbekende u(t,x,y) gegeven op een rand van het integratiegebied.
  • Bij Von Neumann-randvoorwaarden wordt een voorwaarde op het gedrag van de afgeleide van u(t,x,y) op een rand van het integratiegebied gegeven. Meestal is dit de eis dat de partiële afgeleide loodrecht op de rand nul is. Dit stemt fysische overeen met een ondoordringbare rand.

In het geval van meerdere ruimtelijke dimensies dient ervoor gezorgd te worden dat de randvoorwaarden op de hoekpunten onderling niet in tegenspraak met elkaar zijn. Verder moeten ook de randvoorwaarden op het startmoment t = 0 continu aansluiten op de beginvoorwaarden.

Tijd- en ruimtelijke stap[bewerken]

De tijdstap en de ruimtelijke stap kunnen niet los van elkaar gekozen worden. Teneinde de stabiliteit te behouden dient de tijdstap kleiner gekozen te worden dan een bepaalde grens bepaald door de ruimtelijke stap en de diffusiecoëfficiënt k:

dt \, < \, \frac{h^2}{2k} \!

Voorbeeld met één ruimtelijke dimensie[bewerken]

Oplossing van het 1-dimensionale probleem beschreven in het artikel. De effecten van de twee verschillende soorten randvoorwaarden zijn te zien: Rechts (x=0) de Dirichlet voorwaarde die de temperatuur daar op elk moment op nul houdt, links (x=100) de Von Neumannvoorwaarde die een ondoordringbare wand vertegenwoordigt door middel van een ruimtelijke afgeleide die nul is

Probleemstelling[bewerken]

De diffusievergelijking in één ruimtelijke dimensie beschrijft bijvoorbeeld de temperatuur in een dunne staaf met lengte L:

 \frac{\partial u}{\partial t} = k \,  \frac{\partial^2 u}{\partial x^2} \!

De beginvoorwaarde is dan bijvoorbeeld:

u(0,x) \, =  1 \, - \, cos(2\pi x /L)

Dit is een gewone sinusvormige golf die start in x = 0, zijn maximum bereikt in x = L/2 en terug nul wordt x = L.

Als randvoorvaarden kan men bijvoorbeeld op x = 0 als Dirichlet-randvoorwaarde vragen dat de temperatuur daar op elk moment nul is, en op x = L als Von Neumann-randvoorwaarde dat de warmteflux doorheen het grenspunt x = L nul is. Dus respectievelijk:

u(t,0) \, = \, 0 \!
\frac{du}{dx}(t,L) \, = \, 0 \!

Dit alles beschrijft het probleem nu op een unieke manier.

Oplossing door middel van de Lijnmethode[bewerken]

Om de ruimtelijke coördinaat discreet te maken dient de lengte van het lijnstuk [0,L] te worden opgedeeld in N intervallen met lengte h = L/N. Zo ontstaan de netpunten:

x_i \, = \, i.h \!

waarbij

x_0 \, = \, 0 \!
x_N \, = \, L \!

De diffusievergelijking wordt dus, door gebruik te maken van de meest eenvoudige eindige differenties voor de tweede afgeleide:

 \frac{\partial u_i}{\partial t} \, = \, \frac{k}{h^2} \, (u_{i-1} \,  -  2u_i \, +  u_{i+1} \, )

Deze betrekking is geldig voor de index i = 1 .. N-1, terwijl de waarden i = 0 en i = N worden gegeven door de randvoorwaarden. De partiële differentiaalvergelijking is dus omgezet in N-1 gekoppelde gewone differentiaalvergelijkingen. Elke differentiaalvergelijking is gekoppeld aan de differentiaalvergelijkingen horend bij de twee naburige netpunten, of aan een randvoorwaarde en de differentiaalvergelijking van het naburig netpunt.

In het bijzonder geldt voor de beginvoorwaarden:

u(0,x_i) = sin(\pi x_i/L) \!

En voor de randvoorwaarden:

  • de waarde u(0) is overal nul wegens de daar opgelegde Dirichlet-voorwaarde. In dit geval dus:
u_0 \, = \, 0 \!
  • de waarde u(N) wordt bepaald door de Von Neumaanvoorwaarde en kan afgeleid worden uit de achterwaartse tweede orde eindigedifferentie van de eerste afgeleide. Op elk moment t kan de afgeleide in het punt  x_N = L worden benaderd door middel van de functiewaarden in dat punt en in de twee vorige netpunten met index N-1 en N-2. Door deze eindige differentie gelijk aan nul te stellen vindt men:
u_N \, = \, \frac{4u_{N-1}-u_{N-2}}{3} \!

Een andere veelgebruikte en elegante oplossing om een Von Neumannvoorwaarde waarbij de afgeleide op de rand nul is, in te bouwen, is het toevoegen van een extra netpunt N+1, dus buiten het te berekenen interval. De functiewaarde in dat extra netpunt neemt men dan gelijk aan de functiewaarde van het voorlaatste netpunt:

u_{N+1} \, = \, u_{N-1} \!

Doordat in de ruimtelijke richting met een driepunts eindige differentie wordt gewerkt betekent dit dat de afgeleide in het randpunt N automatisch nul wordt. Een parabool door drie equidistente punten, waarvan de twee buitenste dezelfde functiewaarde hebben heeft immers zijn top in het middenste punt, en in de top van de parabool is de afgeleide nul. Dit principe van een extra netpunt kan ook worden gebruikt indien een niet-nulle flux doorheen de wand gevraagd wordt. In dat geval:

u_{N+1} \, = \, u_{N-1}  \, + \, 2.g.h \!

waarbij g de uitgaande(!) flux doorheen de rand is, en h de ruimtelijke stapgrootte.

De N-1 gekoppelde gewone differentiaalvergelijkingen worden vervolgens numeriek opgelost, en op elk moment t worden de bijhorende waarden in de twee eindpunten van de ruimtelijk interval afgeleid uit de randvoorwaarden. Dit dient ook onderweg te gebeuren tijdens het uitvoeren van de numerieke stappen, wanneer een andere methode dan de meest eenvoudige Eulermethode gebruikt wordt. Zo zal de Runge-Kuttamethode gebruik maken van vier schattingen van de afgeleide onderweg, en dus zullen ook de randwaarden van u(t,x) telkens moeten bepaald worden.

Stelsel van partiële differentiaalvergelijkingen[bewerken]

De methode kan ook worden toegepast op stelsels van gekoppelde partiële differentiaalvergelijkingen. Het aantal gekoppelde gewone differentiaalvergelijkingen neemt dan wel evenredig toe. Een voorbeeld hiervan is te vinden in het artikel over TDR-vergelijkingen. In het model dat daar behandeld wordt komt naast diffusie ook taxis voor.

Referenties[bewerken]

  • W.Hunsdorfer, J.G.Verwer (2007) Numerical Solution of Time-Dependent Advection-Diffusion-Reaction equations, Springer, ISBN 978-3-540-03440-7

Externe links[bewerken]

  • (en) S.Hamdi, W.E.Schiesser, G.W.Griffiths, (2007) Method of Lines, Scholarpedia, zie hier
  • (en) W.E.Schiesser, (2000) Partial Differential Equations, [1]