Einführung in Mathematica

Carsten Rezny

Institut für Angewandte Mathematik, Universität Bonn

Eine andere Möglichkeit, eine Funktion in zwei Variablen darzustellen, ist DensityPlot:

In[1]:= ManipulateADensityPlotASinAx2 + y2E, 8x, -a, a<, 8y, -a, a<E, 88a, 1<, 0.01, 5<E

Out[1]=

a

-1.0 -0.5 0.0 0.5 1.0

-1.0

-0.5

0.0

0.5

1.0

Eine andere Möglichkeit, eine Funktion in zwei Variablen darzustellen, ist ContourPlot:

In[2]:= ContourPlotASinAx2 + y2E, 8x, -3, 3<, 8y, -3, 3<E

Out[2]=

-3 -2 -1 0 1 2 3

-3

-2

-1

0

1

2

3

Auch hier gibt es zahlreiche Optionen, um die Details der Ausgabe zu beeinflussen.

In[3]:= ContourPlotASinAx2 - y2E, 8x, -3, 3<, 8y, -3, 3<, Contours ® 4E

Out[3]=

-3 -2 -1 0 1 2 3

-3

-2

-1

0

1

2

3

2 vl4.nb

In[4]:= ContourPlotASinAx2 - y2E, 8x, -3, 3<, 8y, -3, 3<,

Contours ® 30, ColorFunction ® "StarryNightColors"E

Out[4]=

-3 -2 -1 0 1 2 3

-3

-2

-1

0

1

2

3

Es gibt zahlreiche Farbschemata:

In[5]:= ColorData@"Gradients"D

Out[5]= 8AlpineColors, Aquamarine, ArmyColors, AtlanticColors, AuroraColors, AvocadoColors,

BeachColors, BlueGreenYellow, BrassTones, BrightBands, BrownCyanTones,

CandyColors, CherryTones, CMYKColors, CoffeeTones, DarkBands, DarkRainbow,

DarkTerrain, DeepSeaColors, FallColors, FruitPunchColors, FuchsiaTones, GrayTones,

GrayYellowTones, GreenBrownTerrain, GreenPinkTones, IslandColors, LakeColors,

LightTemperatureMap, LightTerrain, MintColors, NeonColors, Pastel, PearlColors,

PigeonTones, PlumColors, Rainbow, RedBlueTones, RedGreenSplit, RoseColors, RustTones,

SandyTerrain, SiennaTones, SolarColors, SouthwestColors, StarryNightColors,

SunsetColors, TemperatureMap, ThermometerColors, ValentineTones, WatermelonColors<Ähnlich wie ContourPlot zeigt auch DensityPlot den Graphen der Funktion in der Draufsicht, verwendet aber statt Höhenlinien einen

kontinuierlichen Farbverlauf.

vl4.nb 3

In[6]:= DensityPlotASinAx2 - y2E, 8x, -3, 3<, 8y, -3, 3<, ColorFunction ® "LightTerrain"E

Out[6]=

-3 -2 -1 0 1 2 3

-3

-2

-1

0

1

2

3

Die Auflösung eines Plots wird durch die Optionen PlotPoints und MaxRecursion bestimmt. Zunächst bestimmt Plot die Funktion-

swerte an PlotPoints Stellen im darzustellenden Bereich. Wenn der Graph zwischen diesen Punkten noch zu “eckig” ist, werden bis zu

MaxRecursion Zwischenwerte bestimmt.

In[7]:= PlotBSin@xD

x

, 8x, -10, 10<, PlotPoints ® 10, MaxRecursion ® 0F

Out[7]=

-10 -5 5 10

0.2

0.4

0.6

0.8

Anschaulich wird das Verhalten mit Manipulate:

4 vl4.nb

In[8]:= ManipulateBPlotBSin@xD

x

, 8x, -10, 10<, PlotPoints ® pp, MaxRecursion ® mrF,

8pp, 2, 100, 1<, 8mr, 0, 10, 1<F

Out[8]=

pp

mr

-10 -5 5 10

-0.10

-0.08

-0.06

-0.04

-0.02

Zurück zur Option PlotRange:

In[9]:= PlotB1

x

, 8x, -5, 5<F

Out[9]=

-4 -2 2 4

-2

-1

1

2

vl4.nb 5

In[10]:= PlotB1

x

, 8x, -5, 5<, PlotRange ® AllF

Out[10]=

-4 -2 2 4

500

1000

Wenn ein Plot beliebig große (oder kleine) Werte enthalten kann, bestimmt die Auflösung, wie weit All geht.

In[11]:= PlotB1

x

, 8x, -5, 5<, PlotRange ® All, PlotPoints ® 200F

Out[11]=

-4 -2 2 4

-1000

1000

2000

3000

4000

5000

Die Funktion ParametricPlot zeichnet parametrisierte Kurven.

Z.B. Kreise

6 vl4.nb

In[12]:= ParametricPlot@8Sin@tD, Cos@tD<, 8t, 0, 2 Π<D

Out[12]=

-1.0 -0.5 0.5 1.0

-1.0

-0.5

0.5

1.0

Ellipsen

In[13]:= ParametricPlot@8Sin@tD, .6 Cos@tD<, 8t, 0, 2 Π<D

Out[13]=

-1.0 -0.5 0.5 1.0

-0.6

-0.4

-0.2

0.2

0.4

0.6

Spiralen

vl4.nb 7

In[14]:= ParametricPlotA9t Log@tD2Sin@tD, t Log@tD2

Cos@tD=, 8t, 0, 8 Π<E

Out[14]=

-200 -100 100

-200

-100

100

200

oder Lissajous-Figuren

In[15]:= ParametricPlot@8Sin@3 tD, Sin@12 tD<, 8t, 0, 2 Π<D

Out[15]=

-1.0 -0.5 0.5 1.0

-1.0

-0.5

0.5

1.0

8 vl4.nb

In[16]:= Manipulate@ParametricPlot@8Sin@3 t + uD, Sin@12 tD<, 8t, 0, 2 Π<D, 8u, 0, 2 Π<D

Out[16]=

u

-1.0 -0.5 0.5 1.0

-1.0

-0.5

0.5

1.0

Zufallszahlen und Statistik

Pseudozufallszahlen können einzeln oder als Liste aus einem frei wählbaren Bereich erzeugt werden.

In[17]:= RandomInteger@100D

Out[17]= 54

In[18]:= RandomInteger@81000, 1999<D

Out[18]= 1091

In[19]:= RandomInteger@81000, 1999<, 10D

Out[19]= 81049, 1452, 1011, 1546, 1626, 1595, 1771, 1976, 1972, 1706<Für Gleitkommazahlen gibt es RandomReal; RandomComplex für komplexe Zufallszahlen. RandomReal arbeitet per default mit einer

Präzision von sechs Stellen, was durch die Option WorkingPrecision beeinflusst werden kann.

In[20]:= RandomReal@10D

Out[20]= 8.06614

In[21]:= RandomReal@10, WorkingPrecision ® 25D

Out[21]= 6.997275582593105916826172

In[22]:= RandomComplex@10 + 10 ä, WorkingPrecision ® 25D

Out[22]= 0.9724387591783456736132996 + 2.270163456404061137116325 ä

Der Zufallszahlengenerator berechnet eine Sequenz von annähernd gleich verteilten Zahlen, die durch ihren Startwert eindeutig bestimmt ist

(daher Pseudo-Zufallszahlen). Der Startwert wird beim Start von Mathematica aus einer möglichst zufälligen Quelle gesetzt (/dev/random,

Systemuhr). Zum Testen einer Rechnung ist es oft nützlich, reproduzierbare Zufallszahlen zu erhalten. Dafür kann man den Startwert manuell

setzen:

vl4.nb 9

Der Zufallszahlengenerator berechnet eine Sequenz von annähernd gleich verteilten Zahlen, die durch ihren Startwert eindeutig bestimmt ist

(daher Pseudo-Zufallszahlen). Der Startwert wird beim Start von Mathematica aus einer möglichst zufälligen Quelle gesetzt (/dev/random,

Systemuhr). Zum Testen einer Rechnung ist es oft nützlich, reproduzierbare Zufallszahlen zu erhalten. Dafür kann man den Startwert manuell

setzen:

In[23]:= SeedRandom@42D; RandomReal@D

Out[23]= 0.425905

In[24]:= RandomReal@D

Out[24]= 0.391023

In[25]:= SeedRandom@42D; RandomReal@D

Out[25]= 0.425905

Die Funktion BlockRandom lokalisiert den Zufallszahlengenerator, so dass die Verwendung von Random-Funktionen innerhalb von

BlockRandom sich nicht auf später generierte Zufallszahlen auswirkt.

In[26]:= 8RandomReal@D, BlockRandom@SeedRandom@42D; RandomReal@DD, RandomReal@D<

Out[26]= 80.391023, 0.425905, 0.347069<Einige Statistikfunktionen:

In[27]:= list = RandomReal@50, 1000D;Mean@listD

Out[28]= 25.1434

In[29]:= StandardDeviation@listD

Out[29]= 14.5994

Anstelle von gleichverteilten Zufallszahlen kann man auch andere Verteilungen verwenden.

In[30]:= RandomInteger@BinomialDistribution@100, 0.7D, 30D

Out[30]= 871, 68, 69, 75, 71, 60, 64, 68, 70, 62, 68, 74, 69,

75, 71, 64, 69, 68, 68, 71, 70, 70, 69, 75, 71, 69, 67, 70, 63, 68<Werte zählen:

In[31]:= % �� Tally

Out[31]= 8871, 5<, 868, 6<, 869, 5<, 875, 3<, 860, 1<,864, 2<, 870, 4<, 862, 1<, 874, 1<, 867, 1<, 863, 1<<

In[32]:= list = RandomInteger@BinomialDistribution@100, 0.3D, 100000D;

In[33]:= ListPlot@Tally@listD, Filling ® AxisD

Out[33]=

20 30 40 50

2000

4000

6000

8000

10 vl4.nb

Oder einfacher:

In[34]:= Histogram@listD

Out[34]=

20 30 40 50

2000

4000

6000

8000

Es gibt eine größere Auswahl an Wahrscheinlichkeitsverteilungen: NormalDistribution, CauchyDistribution,

PoissonDistribution, ...

Die Wahrscheinlichkeitsdichtefunktion (probability density function, PDF) einer Verteilung kann man auch direkt plotten:

In[35]:= DiscretePlot@PDF@PoissonDistribution@15D, xD, 8x, 0, 50<D

Out[35]=

10 20 30 40 50

0.02

0.04

0.06

0.08

0.10

Ein Beispiel für eine einfache Irrfahrt (Random Walk):

In[36]:= data = Module@8t = 0<, Table@t += RandomReal@8-1, 1<D, 82000<DD;

In[37]:= ListPlot@dataD

Out[37]=

500 1000 1500 2000

-10

-5

5

10

15

vl4.nb 11

Lösen von Gleichungen

Viele Gleichungen kann man mit Solve lösen.

In[38]:= eqn = x4

- 4 x3

+ 3 x - 1 � 0

Out[38]= -1 + 3 x - 4 x3

+ x4

� 0

In[39]:= sol = Solve@eqn, xD

Out[39]= ::x ®1

2

2 - 5 - 7 - 2 5 >, :x ®1

2

2 - 5 + 7 - 2 5 >,

:x ®1

2

2 + 5 - 7 + 2 5 >, :x ®1

2

2 + 5 + 7 + 2 5 >>Das Resultat ist eine geschachtelte Liste von Ersetzungsregeln. Dazu ein kurzer

Exkurs: Ersetzungsregeln

Ein Ausdruck der Form x ® y2 ist eine Ersetzungsregel. Im Unterschied zur Zuweisung wird eine Ersetzungsregel nur auf Anforderung

angewandt, nie automatisch.

Die Anwendung einer Ersetzungsregel geschieht mit dem Operator /. (ReplaceAll).

In[40]:= x2

+ 2 x - 1 �. x ® y2

Out[40]= -1 + 2 y2

+ y4

Der Pfeil-Operator ® kam schon bei Funktionsoptionen vor, die tatsächlich auch Ersetzungsregeln sind.

Ersetzungsregeln können auch Muster verwenden:

In[41]:= f@x^2D + f@yD - 2 fAHu + 1L2E �. fAx_2E ® g@xD

Out[41]= f@yD - 2 g@1 + uD + g@xDAnalog zu verzögerten Zuweisungen (:=) gibt es auch bei Ersetzungsregeln eine Variante mit verzögerter Auswertung der rechten Seite (:>).

In[42]:= a = 2; rule = x ® 2 a + 1

Out[42]= x ® 5

In[43]:= rule = x ¦ 2 a + 1

Out[43]= x ¦ 2 a + 1

Ersetzungen werden nicht rekursiv durchgeführt. Wird ./ mit einer Liste von Ersetzungsregeln verwendet, werden die Regeln der Reihe nach

probiert, bis eine passt.

In[44]:= x �. 8x -> 1, x -> 3, x -> 7<

Out[44]= 1

Eine geschachtelte Liste von Regeln erzeugt eine Liste von Resultaten:

In[45]:= x �. 88x -> 1<, 8x -> 3<, 8x -> 7<<

Out[45]= 81, 3, 7<Exkurs Ende

Deshalb ist auch das Resultat von Solve eine Liste von Listen. Wendet man es auf die gelöste Gleichung an, erhält man

12 vl4.nb

In[46]:= eqn �. sol

Out[46]= :-1 +3

2

2 - 5 - 7 - 2 5 -1

2

2 - 5 - 7 - 2 5

3

+1

16

2 - 5 - 7 - 2 5

4

� 0,

-1 +3

2

2 - 5 + 7 - 2 5 -1

2

2 - 5 + 7 - 2 5

3

+1

16

2 - 5 + 7 - 2 5

4

� 0,

-1 +3

2

2 + 5 - 7 + 2 5 -1

2

2 + 5 - 7 + 2 5

3

+1

16

2 + 5 - 7 + 2 5

4

� 0,

-1 +3

2

2 + 5 + 7 + 2 5 -1

2

2 + 5 + 7 + 2 5

3

+1

16

2 + 5 + 7 + 2 5

4

� 0>In[47]:= Simplify@%D

Out[47]= 8True, True, True, True<Für Gleichungen fünften und höheren Grades können die Lösungen im Allgemeinen nicht in geschlossener Form notiert werden.

In[48]:= SolveAx5 + 4 x3

- 1 � 0, xE

Out[48]= 99x ® RootA-1 + 4 ð13

+ ð15&, 1E=,

9x ® RootA-1 + 4 ð13

+ ð15&, 2E=, 9x ® RootA-1 + 4 ð1

3+ ð1

5&, 3E=,

9x ® RootA-1 + 4 ð13

+ ð15&, 4E=, 9x ® RootA-1 + 4 ð1

3+ ð1

5&, 5E==

Dennoch sind diese Lösungen exakte Zahlen, die Mathematica durch sogenannte Root-Objekte darstellt. Mit Root-Objekten kann man wie mit

anderen Zahlen weiterrechnen.

In[49]:= x5

+ 4 x3

- 1 �. x ® RootA-1 + 4 ð13

+ ð15&, 1E

Out[49]= -1 + 4 RootA-1 + 4 ð13

+ ð15&, 1E3

+ RootA-1 + 4 ð13

+ ð15&, 1E5

In[50]:= Simplify@%D

Out[50]= 0

Wenn eine numerische Lösung gesucht ist, kann auch das deutlich schnellere NSolve verwendet werden:

In[51]:= NSolveAx5 + 4 x3

- 1 � 0, xE

Out[51]= 88x ® -0.336836 - 0.542772 ä<, 8x ® -0.336836 + 0.542772 ä<,8x ® 0.0310993 - 2.00169 ä<, 8x ® 0.0310993 + 2.00169 ä<, 8x ® 0.611473<<

In[52]:= Clear@a, b, cD

Generische und nicht-generische Lösungen

Für parametrisierte Gleichungen wie

In[53]:= eqn = a x2

+ b x + c � 0;

liefert Solve nur generische Lösungen, also solche, die sich explizit auf die gesuchte Variable beziehen. Dabei macht Solve implizit

Annahmen über die Parameter.

In[54]:= Solve@eqn, xD

Out[54]= ::x ®-b - b

2 - 4 a c

2 a

>, :x ®-b + b

2 - 4 a c

2 a

>>In diesem Fall liefert Reduce eine allgemeinere Lösung, die auch die vorhandenen Parameter berücksichtigt:

vl4.nb 13

In[55]:= Reduce@eqn, xD

Out[55]= a ¹ 0 && x �-b - b

2 - 4 a c

2 a

ÈÈ x �-b + b

2 - 4 a c

2 a

ÈÈ

a � 0 && b ¹ 0 && x � -c

b

ÈÈ Hc � 0 && b � 0 && a � 0LReduce liefert auch Lösungsbedingungen für Ungleichungen.

In[56]:= ReduceAx3 - 2 x2

+ 1 > 0, xE

Out[56]=

1

2

J1 - 5 N < x < 1 ÈÈ x >1

2

J1 + 5 NIm Unterschied zu Solve gibt Reduce keine Ersetzungsregeln zurück, sondern logische Ausdrücke von (Un-)Gleichungen. Im Fall von

Gleichungen kann man diese mit ToRules wieder in Ersetzungsregeln umwandeln:

In[57]:= ToRules@Reduce@eqn, xDD

Out[57]= SequenceB:x ®-b - b

2 - 4 a c

2 a

>, :x ®-b + b

2 - 4 a c

2 a

>, :a ® 0, x ® -c

b

>, 8c ® 0, b ® 0, a ® 0<FNumerische Lösungen

Wenn weder Solve noch Reduce eine Lösung finden, kann mit FindRoot numerisch nach Lösungen gesucht werden.

In[58]:= eqn = LogAx2 + 1E � x

Out[58]= LogA1 + x2E � x

In[59]:= Reduce@eqn, xD

Reduce::nsmet : This system cannot be solved with the methods available to Reduce. �

Out[59]= ReduceBLogA1 + x2E � x , xF

(Solve und Reduce benötigen unter Umständen recht lange.)

FindRoot sucht ausgehend vom angegebenen Startwert eine Lösung der Gleichung mit dem Newton-Verfahren. Wie viele Lösungen es gibt,

erfährt man so nicht.

In[60]:= FindRoot@eqn, 8x, 0<D

Out[60]= 8x ® 0.<In[61]:= FindRoot@eqn, 8x, 1<D

Out[61]= 8x ® 1.59038<

14 vl4.nb