Institut für Theoretische Physik Universität zu Köln
\n",
"
\n",
"
Prof. Dr. Simon Trebst Christoph Berke
\n",
"
\n",
"
\n",
"\n",
"
Statistische Physik
\n",
"
Übungsblatt 3
\n",
"\n",
"
Wintersemester 2020/21
\n",
"\n",
"\n",
"**Website** [http://www.thp.uni-koeln.de/trebst/Lectures/2020-StatPhys.shtml](http://www.thp.uni-koeln.de/trebst/Lectures/2020-StatPhys.shtml)\n",
"\n",
"**Abgabe**: Montag,23.11.2020, bis 10:00 Uhr\n",
"\n",
"**Besprechung**: Dienstag, 24.11.2020\n",
"\n",
"**Name**: Bitte geben Sie ihren Namen an\n",
"\n",
"**Matrikelnummer**: Bitte geben Sie ihre Matrikelnummer an"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Satz von Liouville\n",
"\n",
"In der folgenden Aufgabe sollen Sie die Phasenraumdynamik verschiedener physikalischer Systeme untersuchen. Das Lösen der kanonischen Gleichungen wird von dem bereits in der Computerphysik benutzten Paket [``DifferentialEquations``](https://github.com/JuliaDiffEq/DifferentialEquations.jl) übernommen. Ihre Aufgabe ist es\n",
"- an den gekennzeichneten Stellen im Programm die Phasenraumvolumina, deren zeitlicher Verlauf untersucht wird, auszuwählen und\n",
"- die Bewegungsgleichungen für den gedämpften harmonischen Oszillator und die auf dem Aufgabenblatt gegebene Hamilton-Funktion zu ergänzen. \n",
"\n",
"Als Hilfestellung ist der harmonische Oszillator bereits implementiert."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"using DifferentialEquations\n",
"using PyPlot\n",
"pygui(true)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"function handle_close(event)\n",
" global update_figure\n",
" update_figure = false\n",
" print(\"Figure closed\")\n",
"end"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Definition der Bewegungsgleichungen für die drei zu untersuchenden Systeme."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"function oscillator(du, u, params, t)\n",
" du[1] = u[2]\n",
" du[2] = -u[1]\n",
"end\n",
"\n",
"function hamiltonian(du, u, params, t)\n",
" # Hier die Bewegungsgleichungen der Hamilton-Funktion auf dem Aufgabenblatt implementieren\n",
"end\n",
"\n",
"function damped_oscillator(du, u, params, t)\n",
" # Hier die Bewegungsgleichungen des gedämpften harmonischen Oszillators implementieren\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"h = oscillator # Wahl des Modells\n",
"\n",
"xs = Array{Float64}[] # Array um die Zeitreihen der Orte zu speichern\n",
"ps = Array{Float64}[] # Array um die Zeitreihen der Impulse zu speichern\n",
"ts = (0. , 50.) # Start und Endzeit der Simulation\n",
"\n",
"\n",
"# Lege Phasenraumvolumen fest und löse die DGL für alle Punkte in dem Volumen.\n",
"\n",
"# Veraendern Sie die Schleife so, dass alle Koordinatenpaare (x,p) in einem Kreis\n",
"# mit Radius 5 um den Ursprung liegen. Als Abstand zwischen den Punkten\n",
"# koennen Sie auf jeder Achse 0.2 auswaehlen\n",
"for x in -1:0.1:0\n",
" for p in 2:0.1:3\n",
" \n",
" init = [x, p] # Anfangsbedingung\n",
" prob = ODEProblem(h, init, ts) # DGL Problem definieren \n",
" sol = solve(prob, RK4(), dt=0.1, adaptive=false) # Problem mit Runge-Kutta 4. Ordnung lösen\n",
" push!(xs, sol[1, :])\n",
" push!(ps, sol[2, :])\n",
" end\n",
"end\n",
"\n",
"# Für eine besonders anschauliche Visualisierung betrachten wir den Rand des Phasenraumvolumens gesondert\n",
"bound_xs = Array{Float64}[]\n",
"bound_ps = Array{Float64}[]\n",
"\n",
"# Veraendern Sie die Schleife so, dass alle Koordinaten in einem Kreis\n",
"# mit Radius 5 < r < 5.2 um den Ursprung liegen. Als Abstand zwischen den Punkten\n",
"# koennen Sie auf jeder Achse 0.05 auswaehlen\n",
"for x in -1:0.1:0\n",
" for p in -1:0.1:0\n",
" init = [x, p]\n",
" prob = ODEProblem(h, init, ts)\n",
" sol = solve(prob, RK4(), dt=0.1, adaptive=false)\n",
" xt = sol[1, :]\n",
" pt = sol[2, :]\n",
" push!(bound_xs, xt)\n",
" push!(bound_ps, pt)\n",
" end\n",
"end"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Animation"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Initiales Phasenraumvolumen\n",
"px = [x[1] for x in xs]\n",
"pp = [p[1] for p in ps]\n",
"\n",
"# Initialer Rand des Phasenraumvolumens\n",
"bpx = [x[1] for x in bound_xs]\n",
"bpp = [p[1] for p in bound_ps]\n",
"\n",
"# Zwei Phasenraumpunkte werden gesondert verfolgt.\n",
"i1 = 1\n",
"i2 = 100\n",
"\n",
"# Fenster erzeugen\n",
"fig = figure(figsize=(5,5))\n",
"fig.canvas.mpl_connect(\"close_event\", handle_close)\n",
"\n",
"# Ausgangssituation plotten\n",
"s, = plot(px, pp, \"o\", alpha=0.2)\n",
"b, = plot(bpx, bpp, \"o\", color=\"black\", markersize=2)\n",
"s1, = plot(px[i1], pp[i1], \"s\")\n",
"s2, = plot(px[i2], pp[i2], \"s\")\n",
"\n",
"axis(\"off\")\n",
"xlim([-10, 10])\n",
"ylim([-10, 10])\n",
"\n",
"update_figure = true\n",
"\n",
"# Schleife fuer Animation\n",
"for i in 1:500\n",
" px = [x[i] for x in xs]\n",
" pp = [p[i] for p in ps]\n",
" bpx = [x[i] for x in bound_xs]\n",
" bpp = [p[i] for p in bound_ps]\n",
"\n",
" s.set_data(px, pp)\n",
" b.set_data(bpx, bpp)\n",
" s1.set_data(px[i1], pp[i1])\n",
" s2.set_data(px[i2], pp[i2])\n",
" pause(0.0001)\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 1.3.1",
"language": "julia",
"name": "julia-1.3"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.3.1"
}
},
"nbformat": 4,
"nbformat_minor": 2
}