Mathematik online lernen im Mathe-Forum. Nachhilfe online

Startseite » Forum » Randwertproblem 2. Ordnung Loesen mit Matlab

Randwertproblem 2. Ordnung Loesen mit Matlab

Universität / Fachhochschule

Gewöhnliche Differentialgleichungen

Tags: Gewöhnliche Differentialgleichungen, Matlab

Lagrange aktiv_icon

17:11 Uhr, 10.03.2010

Eigentlich habe ich eine ganz einfache ODE 2. Ordnung zu loesen und ich dachte Matlab kann das. Allerdings haben alle Beispiele, die ich im Netz gefunden habe andere Anfangsbedingungen. Ich moechte folgendes Problem loesen:

\frac{d C}{d z^{2}} = F (C)

mit den Randwerten

C (- \infty) = 0.249

und

C (\infty) = 0.751

. Ausserdem habe ich noch C(0) = 0.5 gegeben.

F

ist gegeben als

F (C) = \frac{w (1 - C)}{a} - \frac{w C}{a} + \frac{k T * (l n (C) - l n (1 - C))}{a}

Falls es jemanden interessiert: k ist die Bolzmannkonstante, T eine Temperatur, w = 2.2kT und a ist eine Flaeche. Ich habe folgende Werte verwendet

k = 1.38 \cdot 1 0^{- 23}

T = 295

und

a = 7.8 \cdot 1 0^{- 19}

.

Jetzt habe ich das Problem natuerlich ersteinmal als ein System von ODEs geschrieben:

c_{1} ʹ (z) = c_{2} (z)

c_{2} ʹ (z) = F (c_{2})

c_{2} (- \infty) = 0.249

c_{2} (\infty) = 0.751

und

c_{2} (0) = 0.5

Dann habe ich ein m-File geschrieben, das die ersten beiden Gleichungen enthaelt und dann wollte ich es folgendermassen ausfuehren:

>> zspan=[-1000, 1000];
>> [z, c]=ode45('problem', zspan, c0);

Mein Problem ist jetzt, dass ich nicht weiss wie ich

c 0

definieren kann, da ich kein AWP habe. Waere super wenn ihr mir helfen koenntet.

Vielen Dank

Für alle, die mir helfen möchten (automatisch von OnlineMathe generiert):
"Ich benötige bitte nur das Ergebnis und keinen längeren Lösungsweg."

johannes2010 aktiv_icon

08:53 Uhr, 12.03.2010

Hi,

soweit ich das verstanden habe geht das mit keiner NUMERISCHEN Lsg.
Was ODE grob macht: Du gibst ein AWP vor und er geht step-by-step mit einem numerischen Algorithmus unter Verwendung der DGL vor.
Wenn du dein Problem numerisch lösen möchtest, musst du dein kompletten AWP an einer Stelle schätzen. Anschließend die Berechnung durchführen und checken ob auch die anderen vorgegebenen AWP erreicht wurden. Diesen Weg wird man i.A. nicht wählen.

Du brauchst eine analytische Lsg. des Problems. Da würde ich dir Mathematica (Maple) empfehlen.

Ist das kleine

c

in deiner Gleichung das große

C

johannes2010 aktiv_icon

09:12 Uhr, 12.03.2010

Also auf jeden Fall kein stationärer Wert!

Lagrange aktiv_icon

10:39 Uhr, 12.03.2010

Ersteinmal vielen Dank fuer eure Antworten. Mit maple habe ich es ebenfalls versucht (und mit der symbolischen Toolbox in Matlab). Leider gibt mir beides kein Ergebnis. Ich denke mit

C (- \infty)

und

C (\infty)

ist lediglich gemeint, dass die Werte am Rand

C_{a}

und

C_{b}

sind. Ich haette da ganz mutig

C (- 100) = 0.249, C (100) = 0.751

geschrieben.

Ich weiss uebrigens noch dass das System genau eine Loesung hat.

P.S. oh ja Sorry

c = C

ich aendere das sofort - danke fuer den Hinweis!

johannes2010 aktiv_icon

12:40 Uhr, 12.03.2010

Mal umgeschrieben:

F (C) = a_{1} (1 - 2 C) + a_{2} \ln (\frac{C}{1 - C})

Das ist doch schon garnicht für

C \leq 0

definiert,oder?
Wo hast du denn die DGL her?

Lagrange aktiv_icon

13:17 Uhr, 12.03.2010

Nein "C" selbst steht fuer Konzentration und ist daher im Intervall (0,1). Die Logarithmen sichern gerade dass C immer zwischen 0 und 1 bleibt. C ist eine Funktion in z und z selbst ist definiert auf

(- \infty, \infty)

. Daher sind die Randbedinugen schon ok.

Ich glaube zwar leider nicht, dass es dir viel hilft wenn ich dir sage wo das Problem herkommt, mache ich aber trotzdem gerne:

Vielleicht hast du schon mal etwas von der Cahn-Hilliard Gleichung gehoert. Das ist eine nicht lineare partielle DGL 4. Ordung. Mit ihr kann man z.B. berechnen wie sich eine Fluessigkeit aus 2 verschiedenen Molekuelen entwichelt. Wenn man annimmt, dass Molekuel A von Molekuel A angezogen wird und Molekuel A von Molekuel B abgestossen wird. Die Cahn-Hilliard Gleichung (in einer Dimension) sieht folgender massen aus:

c_{t} = - k^{2} c_{x x x x} + f_{x x} (c)

mit periodischen oder Neumann Randbedingungen (*)

f ist hier die freie Energie des Systems und

k ≪ 1

. In meinem Fall diese freie Energie ist die Stammfunktion der oben angegebenen Funktion

F (C)

. Ein schlauer Kerl (Pego - Front Migration in the Nonlinear Cahn-Hilliard Equation). Sagt nun, dass wenn ich das Problem (in meinem ersten Post) loese, dann bekomme ich das Front-Profil zu (*).

johannes2010 aktiv_icon

14:23 Uhr, 12.03.2010

Ich hab jetzt mal mehrere Varianten probiert. kein Erfolg!
Was sich aber anbieten würde, ist eine Taylor-Reihenentwicklung:

C \in [\frac{1}{4}, \frac{3}{4}]

\ln (\frac{C}{1 - C}) =

siehe Bild

\ln (\frac{C}{1 - C}) \approx 4 \cdot (C - \frac{1}{2}) + 32 \cdot {(C - \frac{1}{2})}^{3} + 1536 \cdot {(C - \frac{1}{2})}^{5} + O_{({(C - \frac{1}{2})}^{7})}

F (C) = - \frac{2 w}{a} (C - \frac{1}{2}) + \frac{k T}{a} \ln (\frac{C}{1 - C})

F (C) \approx (4 \frac{k T}{a} - \frac{2 w}{a}) (C - \frac{1}{2}) + \frac{k T}{a} (32 \cdot {(C - \frac{1}{2})}^{3} + 1536 \cdot {(C - \frac{1}{2})}^{5})

Die Approximation hat dann Ordnung 6!!!
Mit diesem Ansatz kann man zumindest die implizite Form:

\int \int \frac{1}{F (C)} d C d C = \int \int d z d z = \frac{1}{2} z^{2} + z + a_{1}

leicht berechnen.
Für eine explizite Darstellung, würde ich sogar nur den linearen Anteil verwenden.

C_{(z)} = p \cdot \sin ({(\frac{2 (w - 2 k T)}{a})}^{\frac{1}{2}} z) + \frac{1}{2}

p \approx \frac{1}{4}

Lagrange aktiv_icon

17:47 Uhr, 12.03.2010

Also ich habe deine Loesung jetzt mal mit maple gezeichnet und ich muss sagen das sieht sehr sehr gut aus. Deine Loesung geht bei 0 durch 0.5 und sie bewegt sich zwischen 0.75 und 0.25. Viel besser an der Loesung ist noch, dass fuer w = 2*k*T die Loesung konstant 0.5 ist. Konntest du nicht wissen, aber das muss auch so sein :-).

du hast mir wirklich sehr weiter geholfen - vielen vielen Dank. Ich werde dich erwaehnen, wenn ich meinen ersten Nobelpreis gewinne :-)

johannes2010 aktiv_icon

17:50 Uhr, 12.03.2010

Danke. War am Anfang schon sehr verwirrend!
Das Problem ist, dass die Randbedingungen

C (\pm \infty)

nicht zu der Logarithmusfunktion passen.

735967

735192

	Status: nicht eingeloggt	Noch nicht registriert?