Spektral- und Pseudospektralmethoden zur Lösung der Fokker ... Fokker-Planck-Gleichung in Chemie...

Fokker-Planck-Gleichung in Chemie und Physik Numerische Losungsmethoden Sturm-Liouville-Probleme

Spektral- und Pseudospektralmethoden zurLosung der Fokker-Planck und der Schrodinger

Gleichung

Daniel Seibel

Universitat des Saarlandes

5. Juli 2016


Gliederung

1 Fokker-Planck-Gleichung in Chemie und PhysikBrownsche BewegungOrnstein-Uhlenbeck-GleichungRayleigh- und Lorentz-Fokker-Planck-GleichungenSpektrallosung der Rayleigh-Gleichung

2 Numerische LosungsmethodenWiederholungFokker-Planck-GleichungSpektralmethoden mit nichtklassischen BasenPseudospektralmethoden mit nichtklassischen Quadraturen

3 Sturm-Liouville-ProblemeSturm-Liouville-ProblemRotation eines starren KorpersHarmonischer OszillatorZweidimensionale SchrodingergleichungAbschließende Worte


Brownsche Bewegung

Brownsche Bewegung

Betrachte Teilsystem von Teilchen der Masse m, das mitTeilchen eines Hintergrundmediums im Gleichgewicht zurTemperatur Tb interagiert

Skalare Kraft F auf das Brownsche Teilchen ist zufallig, nurdie Reibungskomponente Fs kann als stetig angenommenwerden

Setze F (t) = Fs(t) + Fr (t) mit Zeit t und stochastischerKomponente Fr .


Brownsche Bewegung

Langevin-Gleichung

Mit skalarer Geschwindigkeit v und Reibungskoeffizienten αerhalt man

mdv

dt= −αv(t) + Fr (t)

Problem: Fr ist unbekannt!

Losung: Statt einer nicht deterministischen Gleichungbetrachtet man eine stochastische Differentialgleichung.


Brownsche Bewegung

Ornstein-Uhlenbeck-Gleichung

In diesem Sinne lasst sich aus der Langevin-Gleichung dieOrnstein-Uhlenbeck-Gleichung mit W-Dichte P ableiten

∂P(v , t)

∂t= ν

∂

∂v

[vP(v , t) +

kBTb

m

∂P(v , t)

∂v

],

wobei v(t) Geschwindigkeit zu Zeit t, νv Drift- und ν kBTbm

Diffusionskoeffizient und ν Kollisionsfrequenz sind.



Analytische Losung

Fur t →∞ nimmt P eine stationare Form P0 an

P0(v) =

√m

2πkBTbexp

[−mv2

2kBTb

],

die eine Maxwell-Verteilung ist.

Mit der Anfangsbedingung P(v , 0) = δ(v − v0) lasst sich wiefolgt substituieren

τ =e2νt − 1

ν,

u = veντ ,

P(v , τ) = eντQ(u, τ)



Analytische Losung

Damit wird die Ornstein-Uhlenbeck-Gleichung zu

∂Q(u, τ)

∂τ= D

∂2Q(u, τ)

∂u2

mit Diffusionskoeffizienten D = νkBTbm .

Die Losung kann nun mit Fourier-Transformation bestimmtwerden

Q(u, τ) =1√

4πDτexp

[−(u − u0)2

4Dτ

],

wobei u0 durch die Anfangsbedingung v0 gegeben ist.



Analytische Losung

Rucksubstitution gibt

P(x , t) =1√

π(1− e−2νt)exp

[−(x − x0e

−νt)2

1− e−νt

]

mit reduzierter Geschwindigkeit x =√

mv2

2kBTb.



Spektralmethode

Fur die Spektrallosung substituiert man P(x , t) = e−x2g(x , t)

und erhalt

∂g(x , t)

∂t= ν

[−2x

∂g(x , t)

∂x+∂2g(x , t)

∂x2

]Entwicklung in Hermite-Polynome g(x , t) =

∑∞n=0 cnHn(x)

fuhrt zu

∞∑n=0

Hn(x)dcn(t)

dt= ν

∞∑n=0

cn(t)[−2xH ′n(x) + H ′′n (x)

]



Spektralmethode

Die Identitat −2xH ′n + H ′′n = −2nHn liefert zusammen mitKoeffizientenvergleich

dcn(t)

dt= −2nνcn(t),

also

cn(t) = Hn(x0)e−2nνt

2nn!√π.

Damit ist die Spektraldarstellung

P(x , t) = e−x2∞∑n=0

e−2nνt

2nn!√πHn(x0)Hn(x).


Rayleigh- und Lorentz-Fokker-Planck-Gleichungen


Der Kollisionsoperator fur harte Kugeln in derBoltzmann-Gleichung fur ein Gemisch aus zwei Gasen kannfur die Ubergange

γ = Mm → 0 (Rayleigh-Limit)

γ →∞ (Lorentz-Limit)

durch je eine Fokker-Planck-Gleichung beschrieben werden.



Lorentz-Fokker-Planck-Gleichung

Fur γ →∞ erhalt man die Fokker-Planck-Gleichung

∂P(x , t)

∂t=

1

4

∂

∂x

[(2x2 − 3)P(x , t) +

∂

∂x[xP(x , t)]

],

wobei x =√

mv2

kBTbdie reduzierte Geschwindigkeit ist.



Rayleigh-Fokker-Planck-Gleichung

Fur γ → 0 erhalt man die Fokker-Planck-Gleichung

∂P(y , t)

∂t=

∂

∂y

[(y − 3)P(y , t) +

∂

∂y[yP(y , t)]

]wobei y = mv2

kBTbσ0 die reduzierte Energie mit dem

Wirkungsquerschnitt σ0 ist.


Spektrallosung der Rayleigh-Gleichung

Rayleigh-Gleichung

Setzt man P(y , t) = P0(y)g(y , t) mit P0(y) =√ye−y ,

sodass so erhalt man

∂g(y , t)

∂t= y

∂2g(y , t)

∂y2+

(3

2− y

)∂g(y , t)

∂y.

Die Laguerre-Polynome L12n erfullen

yd2L

12n (y)

dy2+

(3

2− y

)dL

12n (y)

dy= −nL

12n (y).




Die Entwicklung in Laguerre-Polynome zur AnfangsbedingungP(y , 0) = δ(y − y0) liefert

cn =n!

Γ(n + 32 )L

12n (y0)e−nt .

Damit ist die Spektraldarstellung der Losung

P(y , t) =√ye−y

∞∑n=0

n!

Γ(n + 32 )L

12n (y0)L

12n (y)e−nt .


Wiederholung

Differentiationsmatrix

Ausgehend von der Spektraldarstellung des Differential-operators d

dx 1. Ordnung bezuglich einer Orthonormalbasis Pn

dnm =

∫w(x)Pn(x)

dPm(x)

dxdx

erhalt man mit den TransformationsmatrizenTnm =

√wmPn(xm) die pseudospektrale Darstellung

Dij =N−1∑n=0

N−1∑m=0

TindnmTmj =√wiwj

N−1∑m=0

P ′m(xi )Pm(xj).


Fokker-Planck-Gleichung


Fur die Langevin-Gleichung

dv

dt= f (v) + g(v)ζ(t) + η(t)

mit ζ multiplikative und η additive Zufallsvariable und f ,gbekannt, erhalt man die Fokker-Planck-Gleichung

∂P(v , t)

∂t=

∂

∂v

[A(v)P(v , t) +

∂B(v)P(v , t)

∂v

].

Zudem soll die Anfangsbedingung P(v , 0) = δ(v − v0) gelten.




A und B sind zeitunabhangige Drift-undDiffusionskoeffizienten

A(v) = f (v) + βg(v)dg(v)

dv,

B(v) = D + βg(v)2

mit Konstanten β und D.

P nimmt fur t →∞ die stationare Form P0 an

P0(v) =1

B(v)exp

[−∫ v

−∞

A(w)

B(w)dw

].


Spektralmethoden mit nichtklassischen Basen


Man substituiert P(x , t) = P0(x)g(x , t) und erhalt

∂g(x , t)

∂t=

1

P0(x)

∂

∂x

[B(x)P0(x)

∂g(x , t)

∂x

]= −A(x)

∂g(x , t)

∂x+ B(x)

∂2g(x , t)

∂x2.

Ist −L der Operator auf der rechten Seite, so hat man

∂g(x , t)

∂t= −Lg(x , t).



Spektralmethode

L ist selbstadjungiert bezuglich dem Skalarprodukt mitGewichtsfunktion P0, wenn folgendes gilt[

P0(x)B(x)∂g(x , t)

∂x

]∞0

= 0.

Betrachte nun in diesem Fall das Eigenwertproblem

Lψ(x) = λψ(x).



Spektralmethode

Sei {Sn} die Polynombasis orthonormal bezuglich P0. Dannhat L die Spektraldarstellung

L(sp)nm =

∫ ∞0

P0(x)Sn(x)LSm(x)dx .

Partielle Integration unter Voraussetzung der Randbedingungergibt sich die symmetrische Form

L(sp)nm =

∫ ∞0

P0(x)Sn(x)1

P0(x)

d

dx

[P0(x)B(x)

dSm(x)

dx

]dx

= −∫ ∞

0P0(x)B(x)S ′n(x)S ′m(x)dx .


Pseudospektralmethoden mit nichtklassischen Quadraturen

Pseudospektralmethode

Sei {Rn} Orthonormalbasis zur Gewichtsfunktion w mit

Rn(x) =

√P0(x)

w(x)Sn(x).

Eingesetzt in die Spektraldarstellung liefert

L(sp)nm = −

∫ ∞0

w(x)B(x)[R ′n(x) + h(x)Rn(x)

][R ′m(x) + h(x)Rm(x)

]dx .

Dabei ist h gegeben durch

h(x) =w ′

2w(x)− P ′0(x)

2P0(x).




Fur die pseudospektrale Darstellung wahlt man eineQuadraturmethode zur Gewichtsfunktion w

L(sp)nm = −

N∑k=1

wkB(xk)[R ′n(xk) + h(xk)Rn(xk)

][R ′m(xk) + h(xk)Rm(xk)

].

Damit erhalt man die Darstellung

L(ps)ij =

N∑m=1

N∑n=1

TinL(sp)nm Tjm

mit Transformationsmatrizen Tin =√wiRn(xi ).




Die erste Summe ist fur festes k in Rn(xk)

h(xk)N∑

n=1

TinRn(xk) = h(xk)N∑

n=1

√wiRn(xi )Rn(xk)

=h(xk)√wk

δik

und in R ′n(xk)

N∑n=1

TinR′n(xk) =

N∑n=1

√wiRn(xi )

N∑`=1

Dk`

√w`wk

Rn(x`)

=Dki√wk.




Fur die zweite Summe uber m erhalt man ahnliche Ergebnisseund zusammenfassen ergibt die Darstellung

L(ps)ij = −

N∑k=1

B(xk) [Dki + h(xk)δki ] [Dkj + h(xk)δkj ] .

Fur den Fall w = P0 verkurzt sich dies zu

L(ps)ij = −

N∑k=1

B(xk)DkiDkj .



Losung des ursprunglichen Problems

Diagonalisierung von L(sp) gibt die normierten Eigen-funktionen ψn zu den Eigenwerten λn. Die Losung derFokker-Planck-Gleichung ist dann gegeben durch

P(x , t) = P0(x)∞∑n=1

ψn(x0)ψn(x)e−λnt ,

beziehungsweise fur L(ps) mit normierten Eigenfunktionen ψn

zu den Eigenwerten λn

P(xk , tj) ≈ P0(xk)N∑

n=1

ψn(x0)ψn(xk)e−λntj .


Sturm-Liouville-Problem

Definition

Das Sturm-Liouville-Problem ist das Eigenwertproblem

Lψn(x) = λnw(x)ψn(x)

mit einer Gewichtsfunktion w > 0 und demDifferentialoperator L gegeben durch

Lf (x) =d

dx

[p(x)

df (x)

dx

]+ q(x)f (x)

mit Diffusionskoeffizient p und Quellterm q.



Voraussetzungen

Seien von nun an p > 0, dpdx und q reellwertig und stuckweise

stetig. Die Eigenfunktionen ψn seien definiert auf einemIntervall [a, b] ⊂ R und erfullen die Randbedingungen

q(a)ψn(a) + p(a)ψ′n(a) = 0,

q(b)ψn(b) + p(b)ψ′n(b) = 0

Diese Randbedingungen stellen sicher, dass L selbstadjungiertbezuglich dem L2-Skalarprodukt ist.




Die Transformation y =∫ √w(x)

p(x) dx fuhrt zu

−d2φn(y)

dy2+ V (y)φn(y) = λnφn(y)

mit dem Potential V

V (y) =q[x(y)]

w [x(y)]+ m[x(y)]

d2

dy2

1

m[x(y)],

wobei m(x) = [p(x)w(x)]−14 und ψn(x) = m(x)φn[y(x)] gilt.

Orthogonale Polynome erfullen diese Gleichungen furbestimme p und q.


Rotation eines starren Korpers

Das Modell

Die Schrodingergleichung zur Bestimmung derRotationsenergie eines zweiatomigen Molekuls ist

−~2

2I

[1

sin θ

d

dθ

(sin θ

dψ(θ)

dθ

)]= Eψ(θ)

mit der quantisierten Energie E und dem Tragheitsmoment I .

Die Substitution x = cos θ fuhrt zu

Hψ`(x) = − d

dx

[(1− x2)

dψ`(x)

dx

]= λ`ψ`(x)

mit dem dimensionslosen Hamiltonoperator H und demEnergieeigenwert E` = λ`

~2

2I .


Rotation eines starren Korpers

Pseudospektrallosung

Fur die Eigenwerte λ` gilt

λ` = `(`+ 1)

mit Eigenfunktionen ψ`, die gerade die Legendre-Polynome P`sind.

Bezuglich dieser Eigenfunktionen hat H Diagonalgestalt unddie pseudospekrale Losung lautet

H(ps)ij =

N∑k=1

(1− x2k )DkiDkj

mit der Differentiationsmatrix (Dij).


Harmonischer Oszillator


Man betrachte die Schrodingergleichung

Hψn(x) = −1

2ψ′′n(x) +

x2

2ψn(x) = λnψn(x),

Die Eigenwerte λn sind

λn = (n +1

2).

Mit welcher Methode erhalt man die besten Ergebnisse?



Differenzenmethode mit Fourier-Polynomen

Man approximiert den Operator H durch eineDifferenzenmethode und stellt die Eigenfunktionen durchFourierpolynome dar.

Damit erhalt man

Hij =1

2(∆x)2

{π2/3 i = j2(−1)i−j

(i−j)2 i 6= j

}+

x2i

2δij .



Konvergenz der Eigenwerte



Eigenfunktion



Lagrange-Mesh-Methode

Eine Pseudospektralmethode mit Hermite-Polynomen liefertdie Darstellung

Hij =

{(4N − 1− 2x2

i )/12 i = j

(−1)i−j[1/(xi − xj)

2 − 1/4]

i 6= j

}+

x2i

2δij .



Tabelle mit Eigenwerten



Ghost Levels

Es treten falsche Eigenwerte auf, die ghost levels genanntwerden.

Dies wird auf nicht exakte Quadratur zuruckgefuhrt. Furn = m = N ist der Integrand von

Vnm =1

2

∫ ∞−∞

w(x)Hn(x)x2Hm(x)dx

ein Polynom vom Grad 2N + 2, die Quadratur der Ordnung Nist aber nur bis 2N + 1 exakt!




Man substituiert ψn(x) = e−x2

2 φ(x) in der Schrodinger-gleichung und erhalt

1

2φ′′n(x)− xφ′n(x) = nφn(x).

Dies ist gerade eine definierende Gleichung der Hermite-Polynome.

Die Losung kann also direkt angegeben werden

H(ps)ij =

1

2

N∑k=1

DkiDkj .

Eigenwerte sind exakt λn = n und Eigenfunktionenφn(xi ) = Hn(xi ).



Eigenfunktion


Zweidimensionale Schrodingergleichung

Spektralmethode

Der zweidimensionale Hamiltonoperator hat die Form[− ∂2

∂x2− ∂2

∂y2+ V (x , y)

]ψnm(x , y) = λnmψnm(x , y).

Mit Polynombasen Xn(x) zur Gewichtsfunktion u(x) undYm(y) zu v(y) erhalt man analog zum eindimensionalen Fall

Hn′m′,nm = δm′m

∫u(x)X ′n′(x)X ′n(x)dx +

δn′n

∫v(y)Y ′m′(y)Y ′m(y)dy +

(Vn′m′,nm − Vn′m′,nm).



Spektralmethode

Dabei ist V gegeben durch

Vn′m′,nm = δm′m

∫ (1

4U2(x)− 1

2U ′(x)

)u(x)X ′n′(x)X ′n(x)dx +

δn′n

∫ (1

4V 2(y)− 1

2V ′(y)

)v(y)Y ′m′(y)Y ′m(y)dy

mit U(x) = −u′(x)u(x) und V (y) = − v ′(y)

v(y) .

V ist dabei

Vn′m′,nm =

∫ ∫u(x)v(y)Xn′(x)Ym′(y)

V (x , y)Xn(x)Ym(y)dxdy .




Entsprechend ist die pseudospektrale Darstellung

Hij ,k` = δj`

N∑n=1

DniDnk + δik

M∑m=1

DmjDm`

+[V (xi , yj)− V (xi , yj)

]δikδj`

mit

V (x , y) =

(1

4U2(x)− U ′(x)

)+

(1

4V 2(y)− V ′(y)

).



Plot der Eigenfunktionen


Abschließende Worte

Zusammenfassung

Numerisches Losen einer Fokker-Planck-Gleichung kann aufein Eigenwertproblem vom Sturm-Liouville-Typ zuruckgefuhrtwerden.

Spezielle orthogonale Polynome losen das Eigenwertproblem.

Spektral-und Pseudospektrallosungen konnen mit Hilfe dieserorthogonalen Polynome beschrieben werden.


Abschließende Worte

Vielen Dank fur IhreAufmerksamkeit!

Spektral- und Pseudospektralmethoden zur Lösung der Fokker ... Fokker-Planck-Gleichung in Chemie...

Documents

Transcript of Spektral- und Pseudospektralmethoden zur Lösung der Fokker ... Fokker-Planck-Gleichung in Chemie...