Einführung in Mathematica  
-- Teil III  .-
zur Vorlesung Mathematische Methoden der Physik im WiSe 2014/15

: :    Prof. Dr. Norbert Dragon und PD Dr. Michael Flohr    : :    05. 12. 2014    : :

Interaktion

Eine der großen Unterschied von Mathematica zu gewöhnlichen Programmiersprachen ist, dass man sehr einfach interaktive Sandkästen bauen kann, in denen man spielen und ausprobieren kann. Das ist ein unschätzbar wertvolles Hilfsmittel, um eigenständig bei mathematischen Problemen Lösungsstrategien zu finden, oder Ergebnisse in Abhängigkeit von den Parametern zu untersuchen.

Manipulate stellt eine Art "Schachtelfunktion" dar, mit welcher sich ein oder mehrere "innere" Parameter eingebetteter Funktionen (beispielsweise Plot) interaktiv in gegebenen Grenzen variieren lassen. 
Achtung: verwendet man in Manipulate extern definierte Funktionen, dann müssen die Parameter in der externen Definition als Variablen deklariert werden. Sonst behandelt Mathematica diese als Konstanten, und Manipulate kann diese dann nicht mehr variieren.

Manipulation eines Parameters in diskreten Stufen

Manipulate[Factor[x^n-1],{n,1,100,1}]

Manipulation zweier kontinuierlicher Parameter

Manipulate[Plot[Sin[a x+b],{x,0,6}],
{{a,2,"Wellenzahl"},1,4},
{{b,0,"Phase"},0,10}]

Ein komplexeres Beispiel anhand des elektrostatischen Potentiales zweier Ladungen

Manipulate[ContourPlot[q1/Norm[{x,y}-p[[1]]]+q2/Norm[{x,y}-p[[2]]],{x,-2,2},{y,-2,2},Contours->10],
  {{q1,-1},-3,3},
  {{q2,1},-3,3},
  {{p,{{-1,0},{1,0}}},{-1,-1},{1,1},Locator},
  Deployed->True]

Import und Export von Daten

Wir können in Mathematica leicht Daten im- und exportieren. Das kann vor allem dann nützlich sein, wenn man Messdaten mit Mathematica weiterverarbeiten möchte, oder Ergebnisse von Mathematica in spezieller Software für numerische Aufgaben einspeisen möchte (weil Mathematica die symbolischen Ausdrücke nicht mehr weiter lösen kann). Hierbei kann man beim Im- und Export kontrollieren, in welchem Format dies geschehen soll. In unserem ersten Beispiel handelt es sich um Listen, die in Mathematica mit dem Table-Konstrukt erzeugt werden (siehe Teil II).

Export["hurz.txt",Table[{n,N[((((-1))^n)/(n!))]},{n,0,20}],"Table"]

FilePrint["hurz.txt"]

Import["hurz.txt","Table"]

Export["zufall.dat",Table[RandomInteger[{-1000,1000}],{i,1,30},{j,1,30}],"List"];

FilePrint["zufall.dat"]

Statt Import muß nun ReadList benutzt werden, da Mathematica sonst jede Zeile fälschlich als Zeichenkette einliest und dann nicht wie erwartet eine Liste aus dreißig Listen mit jeweils dreißig Elementen erzeugt! Hier wollen wir mal sehen, wie zufällig der Zufallsgenerator von Mathematia würfelt. Die Funktion ListDensityPlot interpoliert die diskreten Werte über einem zwei-dimensionalen Gitter, wobei hellere Farböne größere Werte anzeigen. Sie ist besonders geeignet, die Einträge von Listen von Listen (also einer Liste von Zeilen, so dass in unserem Beispiel ein 30x30 Gitter entsteht) auf einen Schlag zu visualisieren. Je ungeordneter dieses Muster aussieht, desto zufälliger sind die Einträge gewesen.

ListDensityPlot[ReadList["zufall.dat"]]

Export["zufall.jpg",%,"JPEG"]

Anzeigen der zuvor gespeicherten Datei

Import["zufall.jpg"]

Nur zum Vergleich: das folgende Beispiel zeigt ganz klar, dass die Werte in dieser Matrix nicht rein zufällig sein können.

daten = Table[N[Sin[i j]],{i,1,30},{j,1,30}];
ListDensityPlot[daten]

Mit Mathematica sind auch Bildverarbeitungen und -analysen möglich - und zwar mit beliebigen "mathematisch exakten" Filtern. Damit kann man zum Beispiel komplexe Daten sehr gut durch geschickte FIlterung und Visualisierung so darstellen, dass man Zusammenhänge, also mögliche Korrelationen, sichtbar machen kann.

(* Das folgende Kommando müssen Sie auf Ihre Situation anpaasen. 
 Je nach dem, wo Sie die Datei Mainhattan.bmp abspeichern,
   müssen Sie anschließend den Pfad zu Laden des Bildes entsprechend 
  setzen. *)
SetDirectory["mathematica/mmdp/WiSe14-15"]

bild=Import["Mainhattan.bmp"(* by Martin Paech, 2009/02/18 18:46:33 CET *)]

Zunächst müssen die Bilddaten in eine mehrdimensionale Liste umgewandelt werden, ...

bilddaten=ImageData[bild,"Byte"];

bildgroesse=Dimensions[bilddaten]

bildgroesse[[1]]*bildgroesse[[2]]

... welche aus 1167360 Bildpunkten (912 Zeilen zu jeweils 1280 Spalten) besteht, mit einem "Tripel" (dreikomponentiger Vektor) aus ganzen Zahlen für die drei Grundfarben (in 256 Stufen: 0 ist schwarz, 255 entspricht reinem Rot, Grün oder Blau) im RGB-Farbraum.

Um aus dem Negativbild ein Positiv zu machen, ist nur folgendes einfaches Umrechnungsprogramm nötig:

For[y=1,y≤bildgroesse[[1]],y++, (* Achtung: Mathematica vertauscht die "natürlichen" x- und y-Koordinaten! *)
 For[x=1,x≤bildgroesse[[2]],x++,
 bilddaten[[y,x]][[1]]=255-bilddaten[[y,x]][[1]] (* roter Farbkanal *);
 bilddaten[[y,x]][[2]]=255-bilddaten[[y,x]][[2]] (* grüner Farbkanal *);
 bilddaten[[y,x]][[3]]=255-bilddaten[[y,x]][[3]] (* blauer Farbkanal *)
]
]

Man beachte, dass wir diese Manipulation “in palce” vorgenommen haben. Wir erkennen nunmehr, daß Frankfurt am Main im Abendlicht (mindestens ;-) ebenso reizvoll wie Hannover sein kann:

Image[bilddaten,"Byte"]

Doch halt! Dom, Fernsehturm, Großer Feldberg ... befinden sich nördlich des Mains. Auch das ist schnell repariert:

For[y=1,y≤bildgroesse[[1]],y++, 
bilddaten[[y]]=Reverse[bilddaten[[y]]]
]

Image[bilddaten,"Byte"]

Animationen

Nachdem wir gelernt haben, wie wir selbst Daten interaktiv manipulieren können, und dass man sozusagen fast alles zu Daten machen kann, die mit mathematischen und symbolischen Methoden bearbeitet werden können, und nachdem wir vielleicht lange an der Lösung eines Problems gebastelt haben, wollen wir das ganze dann am Ende möglichst cool und effektvoll in Szene setzen - also am besten gleich als eine Art Movie. In der Tat bietet Mathematica auch erstaunlich einfach zu handhabende Möglichkeiten an, Formeln und Daten in animierte, also bewegte, Grafiken umzusetzen.

Clear[x,y]

Das folgende stellt einen brauchbare Approximation für “digitale” Schaltsignale (Pulse) dar. Diese sind selten echte Rechteckfunktionen, sondern sehen durchaus eher so aus (oder noch “schlimmer”) aus.

tp[x_]:=If[(( π)/2)≤ Abs[x]≤((3π)/2),1-Cos[x],0]

Plot[tp[x],{x,-10,10}]

Das folgende Integral nennt man Faltungsintegral. Unser Doppelpuls wird mit einer Gaußkurve gefaltet. Dies simuliert zum Beispiel, wie ein Puls durch Dämpfung in Leitungen seine Form verlieren kann, oder wie zum Beispiel ein eng fokussierter Lichtpuls durch Streuung auseinanderläuft. Ein weiteres typisches Beispiel ist der Verlust von Signalkohärenz durch temperaturbedingte Effekte (thermische Störungen bzw. thermisches Rauschen). 
In unserem Beispiel wollen wir uns vorstellen, dass wir ein Material an zwei Stellen sehr heiß gemacht haben, und wollen sehen, wie es sich mit der Zeit abkühlt. Die Gaußkurve ist nämlich auch hierfür die richtige Funktion, sie spielt bei der Lösung der Wärmeleitgleichung eine wichtige Rolle. Wir haben hier die Breite der Gaußkurve als Funktionsparameter definiert, so dass wir die Auswirkung unterschiedlich starker Streuung bzw. unterschiedlich hoher Temperatur jetzt untersuchen können.

tptd[x_,t_]=(∫_(-∞)^∞)((E^(-((((-x+y))^2)/(4 t))))/(2 (Sqrt[π t])))tp[y]\[DifferentialD]y;

Manipulate[Plot[Re[tptd[x,t]],{x,-10,10},
  AxesLabel->{"x","Temperatur"},
  AxesOrigin->{0,0},
  PlotRange->{0,2}],
{{t,1},(10^(-6)),10}]

Und jetzt wollen wir das ganze als einen Film sehen, wie unsere beiden weißglüheden Stellen im Eisen abkühlen und den restlichen Barren dabei erwärmen (zumindest ist dies eine Interpretation des Vorgangs, die mit der gewählten Visualisierung gut zusammen passt).

Export["Temperaturprofil.avi",Table[DensityPlot[Re[tptd[x,t]/2],{x,-10,10},{y,0,1},AspectRatio->1/8,ColorFunction->"SunsetColors",ColorFunctionScaling->False,Frame->None],{t,(10^(-6)),5,0.1}],"AVI"]

Import["Temperaturprofil.avi","Animation"]

Fouriertransformation

Nach diesen Spielereien wollen wir jetzt zum Abschluss noch ein weiteres wichtiges Anwendungsgebiet von Mathematica vorstellen, nämlich Transformationen von Funktionen. Die wichtigste solche Transformation, die Physiker sehr viel gebrauchen, ist die Fouriertransformation. Hinter ihr steckt eigentlich nicht viel anderes, als die vielleicht intuitiv nachvollziehbare Tatsache, dass ich Musik als Amplitude in der Zeit aufnehmen kann (wie auf der Schallplatte), oder als Frequnzspektrum (wie, etwas vereinfacht betrachtet, im mp3-Spieler). Die Fouriertransformation rechnet die eine Forn in die (völlig äquivalente) andere um.
Auch zwei-dimensional kennen Sie das alle. Um digitale Bilder abzuspeichern, kann man entweder die Helligkeitswerte (Amplituden) eines Farbkanals für die Bildpunkte (zwei-dimensional statt ein-dimensionale Zeit) abspeichern, oder äquivalent die Fouriertransformation davon, was die Frequenzen sind, mit denen sich die Helligekitswerte in verschiedenen Richtungen “wiederholen”. Das ist das, was in einer jpg-Datei passiert, die aus den meisten Digitalkameras rauskommt. Das erstaunliche ist, dass Bilder der Realität nach Fouriertransformation oft erstaunlich gut kompirimiert werden können, so wie auch Musik im mp3-Format bei tolerablen Einbußen in der Klangtreue viel weniger Speicherplatz einnimmt. Übungsaufgabe: Überlegen Sie, woran das liegen könnte ...

Elementare Funktionen

Zunächst ein-dimensional:

l={1,x,(x^2),(1/x),(1/(x^2)),(E^(-((x^2)/2)))};

TableForm[{l,FourierTransform[l,x,k]},TableDirections->Row,TableSpacing->{2,6},TableAlignments->{Center,Left},TableHeadings->{{"f(x)","F(k)"}}]

Und nun zwei-dimensional:

g[x_,y_]:=Exp[-((((x^2)+(y^2)))/2)]

Plot3D[g[x,y],{x,-3,3},{y,-3,3},PlotRange->All,MeshFunctions->{#3&},Mesh->20,ColorFunction->"GreenBrownTerrain"]

FourierTransform[g[x,y],{x,y},{u,v}]

Periodische Funktionen

Strikt periodische Funktionen (also zum Beispiel kontinuierliche Töne) lassen sich durch Fourierreihen darstellen, denn in ihrem Frequenzspektrum treten nur die echten vielfachen der Grundfrequenz auf (Obertöne ...). Wir betrachten nun folgende reihenentwickelte Funktion:

f[t_]:=(∑_(n=0)^7)(1/(2^n))Sin[(2^n)220 2π t]

Plot[f[t],{t,0,0.01}]

Das kann Mathematica auch hörbar machen kann (220 Hz Grundton mit sieben Obertönen, was einen leicht orgelpfeienähnlichen Klang ergibt).

Play[f[t],{t,0,2},SampleRate->22050]

Die (auch aus obiger Formel ersichtlichen) diskreten Frequenzen und Gewichte der Grundschwingung und der Oberschwingungen finden sich als δ-Funktionen in der Fouriertransformierten, der Lösung des Fourierintegrales
	(Sqrt[((|b|)/(((2π))^(1-a)))])(∫_(-∞)^∞)f(t) (e^(i b ω t))dt ,
in welchem die Parameter a und b nachfolgend derart gewählt sind, daß sich unmittelbar das Frequenzspektrum und nicht die "Winkelgeschwindigkeiten" ergeben.

a=0;

b=-2π;

("Hintransformation")

ff[ω_]=FourierTransform[f[t],t,ω,FourierParameters->{a,b}]

("Rücktransformation")

InverseFourierTransform[ff[ω],ω,t,FourierParameters->{a,b}]

Bei der Rücktransformation kann in diese Falle auch auf die Funktion InverseFourierSinTransform zurückgegriffen werden, welche sich auf ungerade, asymmetrische Komponenten beschränkt. Damit vermeiden wir, uns mit der komplexen e-Funktion auseinandersetzen zu müssen.

InverseFourierSinTransform[ff[ω],ω,t,FourierParameters->{a,b}]

Die Nachbarbeitung der allgemeinen Lösung ...

fff[t_]=Simplify[Im[%],Assumptions->t∈Reals]

... führt erwartungsgemäß auf die Ausgangsfunktion ("Beweis durch Hingucken").

Plot[{f[t],fff[t]},{t,0,0.01},PlotStyle->{Directive[Thick,Red],Directive[Dashed,Thick,Black]}]