IWR, Universit¨at Heidelberg Sommersemester 2009
Aufgabenblatt 7 09. Juni 2009
Ubungen zur Vorlesung¨ Simulationswerkzeuge
Dr. S. Lang, D. Popovi´c Besprechung am 16. Juni 2009 in der ¨Ubung
Ubung 13¨ DUNE-pdelab: Einfache Neuronenzellen mit station¨arer Diffusion-Reaktion In der letzten ¨Ubung haben wir die Diffusions-Reaktions-Gleichung in 2D untersucht. Mit dieser Gleichung l¨aßt sich die elektrische Potentialverteilung in einer einfachen zylindrischen Neuronenzelle modellieren. Sei eine idealisierte Neuronenzelle in Form eines Zylinders der L¨angelmit Radiusr, beide inµm, gegeben. Dann beschreibt in der Problemstellung
−πr2GA∂xxv(x) + 2πrGMv(x) = 0 in Ω = [0, l],
∂xv(x) =−2GM
rGA
I0 f¨urx= 0, zus¨atzliche Randbedingungen bei x=l
v(x) die Potentialverteilung bei einem am linken Terminal appliziertem Strom. Weiter bezeichnet GA die spezifische axiale Leitf¨ahigkeit des Zell-Materials (Zytoplasma), w¨ahrend GM die spezifische Leitf¨ahigkeit der Zellmembran ist. Sie haben jeweils die Einheiten Ωcm1 , w¨ahrendI0 inµAangegeben wird.
F¨ur die Randbedingungen am rechten Terminal (x=L) unterscheiden wir zwei F¨alle:
1. Killed end v= 0, entspricht einem Durchtrennen der Zelle,
2. Sealed end ∂xv= 0, entspricht einem Abschluss der Zelle mit einer undurchl¨assigen Membran.
F¨ur die oben angegebenen Probleme existieren analytische L¨osungen. Die allgemeine L¨osungen einer Gleichung der Art a∂xxV −V = 0 l¨asst sich schreiben als v = c1cosh(√xa) +c2sinh(√xa). Mit den Randbedingungen ergeben sich mit 2GrGA
M =:adie L¨osungen:
1. F¨ur den Fall killed end:v(x) = √I0a · n
tanh l cosh
x
√a
−sinh
x
√a
o, 2. F¨ur den Fall sealed end:v(x) = √I0a·
n
coth l cosh
x
√a
−sinh
x
√a
o . Details zu der Herleitung der Gleichung und L¨osung finden sich zum Beispiel in
”Introduction to theoretical neurobiology“, von Henry C. Tuckwell, erschienen beiCambridge University Press.
In Ihrem Home-Verzeichnis ist wie gehabt das Verzeichnis dunemit dem Moduldune-pdelabvor- handen, in dessen Unterverzeichnis dune/pdelab/testwir wieder arbeiten wollen. Zun¨achst basteln wir uns ein neuesmake-Target:
• Dazu f¨ugen wir im Makefile.am ein neues Zielsimw_simpleneuron1D hinzu und kopieren die vier Dateiensimw_diffusionreaction*.{cc,hh}in entsprechende Dateien mit dem neuen Na- men.
• Das Kommando make clean && make simw_simpleneuron1D erzeugt anschließend das Pror- gramm.
Nun m¨ussen wir das Programm problemspezifisch ver¨andern. Dazu gehen wir wie folgt vor:
• Hinzuf¨ugen zweier Klassen f¨ur die exakten L¨osungen der beiden Randwert-Probleme in der Datei simw_simpleneuron1Dfcts.hh. Die Klassen sollen eine Methode T evaluate(X x) ent- halten, die als Parameter eine Koordinatexerh¨alt und einTzur¨uckliefert. Die TypenT,Xsollen Template-Parameter der Funktion sein.
• Inkludieren der im Verzeichnis bereistehenden Dateil2norm.hh. In dieser befindet sich eine tem- platisierte Funktion, die dieL2-Norm des Fehlerse(x) :=v(x)−vh(x) bei gegebener analytischer L¨osungv(x) und approximierter L¨osung vh berechnet.
• Erweiterung der Funktioncylinder1Dum die Ausgabe des Fehlers in derL2-Norm. Dazu wird die Funktionl2normmit entsprechenden Template-Parametern aufgerufen.
Abschließend wollen wir unsere Implementierung testen und dazu folgende zwei Testprobleme l¨osen:
1. Killed end,l= 316µm,r = 2µm,RM = 2.0kΩcm2,GA= 200Ωcm1 ,I0 = 2µA.
2. Sealed end, gleiche Parameter wie bei vorigem Beispiel.
Rechnen Sie beide Probleme und visualisieren Sie die L¨osungen. Das Gleichgewichts-Potential sollte vom linken Terminal weg abfallen. Im Falle killed end ergibt sich nat¨urlich v = 0 am rechten Ter- minal, w¨ahrend der Potentialverlauf bei sealed end Randbedinungen ¨uber dem Potentialverlauf des killed end Falles liegt und am rechten Terminal gr¨oßer 0 ist (Fluß ¨uber den Rand). Kontrollieren Sie die KonvergenzordnungO(h2) des numerischen L¨osungsverfahrens durch Berechnung des Fehlers f¨ur sukzessive verfeinerte Gitter.
