8. Mathematica を使ってみよう
8.1. Mathmatica の特徴
数式処理 数値ではなく，数式を使った演算が可能である． グラフィックス 数式をすぐにグラフ化できる プログラミング 今回の授業では扱わないが，プログラミングが可能である．（プログラミングの文法は C 言語とほぼ同じで す．C 言語をある程度しっていて興味のある人は，試してみてください）8.2. クイックスタート
8.2.1. Mathematica の起動
Mathematica を起動するには，[スタート]メニュ ー か ら [ プ ロ グ ラ ム ] → [Mathematica4] → [Mathematica4]をクリックします．すると右図の ように Mathematica が起動します．8.2.2. 数式の入力と評価
まず，「１＋１」と入力して，[Shift]キーを押しな がら[Enter]キーを押してみましょう．初めだけ 少し時間がかかりますが， In[1]:= 1+1 Out[1]= 2 と表示されるはずです． （注意） ・ 入力した数式の評価を行うためには，常に[Shift]キーを押しながら[Enter]キーを押します．単に[Enter]キ ーを押すと前の行からの入力の続きと判断されます． ・ 「 In[1]:= 」等は Mathematica が勝手に表示するもので入力してはいけません． 次に「x+x」と入力して，[Shift]キーを押しながら[Enter]キーを押して見ましょう．今度は， In[2]:= x+x Out[2]= 2 x と表示されるはずです．このように Mathematica は数値以外の数式を扱うことができます．8.2.3. 計算結果の保存
計算結果を初めて保存する時は，プルダウンメニューの[File]から[Save As…]をクリックします．すると図 8-1 の [名前を付けて保存]ダイアログが開くので，適当な名前をつけて（図 8-1 の例では，「test1」にしてあります）， [保存する場所]で自分のホームディレクトリである[Nafs の Sxxxxx（Y:）]が選択されていることを確かめてから，[保存]ボタンをクリックします． ２回目以降，保存する場合はプルダウンメニューの[File]から[Save]をクリックします． 図 8-1
8.2.4. Mathematica の終了
Mathematica を終了するためには，プルダ ウンメニューの[File]から[Exit]をクリックす るか，タイトルバー（右図）の左上の[X]ボタンをクリックします．8.3. Mathematica のコマンド
8.3.1. 演算コマンド
四則演算 足算，引き算，掛け算，割り算はそれぞれ‘+，-，*，/’で行える．また( )も使うことができる．但し‘[ ]’は関数の 引数とみなされるので使ってはいけない． (例) 1+1，1/3，2*x，x/y，(1+2)*3 冪級数は‘＾’を使う (例) 2^4 その他の組み込み関数 関数 Mathematica の表記 関数 Mathematica の表記x
Sqrt[x]log x
Ｌｏｇ［ｘ］e
x Ｅｘｐ［ｘ］sin x
1 ArcSin[x]sin x
Ｓｉｎ［ｘ］ 1cos x
ＡｒｃＣｏｓ［ｘ］cos x
Ｃｏｓ［ｘ］ 1tan x
ＡｒｃＴａｎ［ｘ］tan x
Ｔａｎ［ｘ］| |
x
Ａｂｓ［ｘ］ Mathematica は大文字・小文字を区別する．組み込み関数は全て大文字で始まることに注意． 組み込まれている定数 定数 Mathematica の表記 定数 Mathematica の表記 Pi ∞ Infinity
8.3.2. 2 次元グラフィックス
陽関数のグラフ Plot[f,{x,xmin,xmax}(,option-> value)] これは，x の関数 f を x の値が xmin から xmax までプロットすることを示します．‘（ ）’内は作図のオプションで 省略できることを示します． (例) ｓｉｎx のプロットPlot[ Sin[x], {x,0,2 Pi} ]
（例） 上と同じで，y 方向の範囲を-2 から２にした例 Plot[ Sin[x], {x,-Pi,Pi},
PlotRange->{-2,2} ]
(例) 外枠をつけて，x，y というラベルをつけた例 Plot[ Sin[x], {x, 0, 2*Pi},
Frame->True,FrameLabel->{x,y} ] -3 -2 -1 1 2 3 -2 -1.5 -1 -0.5 0.5 1 1.5 2 [問題１] (1) x2
e
y
のグラフを
4
x
4
で書きなさい． (2)x
x
y
sin
のグラフを
4
x
4
で書きなさい（グラフ全体が入るようにすること）． 陰関数のグラフ
陰関数をプロットするコマンドは，
ContourPlot[ f, {x, xmin, xmax}, {y,ymin,ymax},(options) ] です． （例）
x
2
y
2
1
のグラフ ContourPlot [x^2+y^2 == 1, {x, -2, 2},{y,-2,2} ] プロットする数式の等号は 2 つ重ねる（「= =」と書く） ことに注意．コンピュータでは「=」は「右辺を左辺に代 入する」という意味で使うのでこれと区別するために， 等号を重ねて書く． -2 -1.5 -1 -0.5 0.5 1 1.5 2 -2 -1.5 -1 -0.5 0.5 1 1.5 2 [問題２] (1)y
2
x
4(
1
x
2)
のグラフを
2
x
2
で書きなさい． (2)x
4
x
2y
2
6
x
2y
y
2
0
のグラフを
4
x
4
で書きなさい． (3)x
4
4
x
2
2
y
2
y
4
1
のグラフを
4
x
4
で書きなさい． (4)x
4
2
x
2
2
y
2
y
4
1
のグラフを
4
x
4
で書きなさい． 2 次元アニメーションの作成 アニメーションを作るためのコマンドは， Animate[変数を含んだ描画関数，{変数名，変数の最小値，変数の最大値，変数の増分}] です． （例） 半径が 0 から１まで 0.1 づつ増えていく円のアニメーションを作る Animate[ContourPlot[x^2+y^2==a^2,{x,-2,2},{y,-2,2}], {a,0,1,0.1}] 解説: ContourPlot[x^2+y^2==a^2,{x,-2,2},{y,-2,2}] この部分はすでに勉強した円を各部分です．ただし半径は a で変数に なっています． Animate 関数の中で半径を 0 から 1 まで 0.1 づつ増や して行きます．それを指定しているのが， {a,0,1,0.1} の部分です． このままでは，図がカクカクするのでもう少し描画するメッシュを細かくし ましょう．そのためのオプションが「PlotPoints」です．Animate[ContourPlot[x^2 + y^2 = = a^2, {x, -2, 2}, {y, -2, 2}, PlotPoints -> 100], {a, 0, 1, 0.1}] -2 -1.5 -1 -0.5 0.5 1 1.5 2 -2 -1.5 -1 -0.5 0.5 1 1.5 2 -2 -1.5 -1 -0.5 0.5 1 1.5 2 -2 -1.5 -1 -0.5 0.5 1 1.5 2
[問題３] 陽関数のアニメーション．波の伝播(
f
(
x
t
)
でかける関数は波の伝播を示すことを目で確かめる) (1)sin
(
x
t
)
のグラフを
10
x
10
，
1
.
5
y
1
.
5
の範囲で時間 t を 0 から 10 まで 0.5 づつ増やし た時の変化をアニメーションにしなさい． (2) (x t)2e
のグラフを
2
x
2
，
1
.
5
y
1
.
5
の範囲で時間 t を 0 から 4 まで 0.2 づつ増やした時の変 化をアニメーションにしなさい． [問題 4] 陰関数のアニメーション （1） 2 3 2ax
x
y
のグラフを
4
x
4
，
5
y
5
の範囲で a を-2 から 2 まで 0.1 づつ増やした時の変 化をアニメーションにしなさい． （2） 2 2 2 2 24
)
1
(
x
y
x
b
のグラフを
2
x
2
，
2
y
2
の範囲で b を 0.1 から 2 まで 0.1 づつ 増やした時の変化をアニメーションにしなさい．8.3.3. 3 次元グラフ
3 次元グラフ
Plot3D[f,{x,xmin,xmax},{y,ymin,ymax},(,option-> value)]
これは，x，ｙの関数 f を x の値が xmin から xmax まで， ｙの値が ymin から ymax までの範囲で計算し z=f(x,y) のグラフをかきます．‘（ ）’内は作図のオプションで省略できることを示します． (例) z=sinx*cosx を-< x <，-< y < の範囲でプロットす る． Plot3D[Sin[x]*Cos[y],{x,-Pi,Pi},{y,-Pi,Pi}] 密度グラフ DensityPlot[f,{x,xmin,xmax},{y,ymin,ymax}(,option->value)] 密度グラフは値の大きさを色の濃淡で表すことで 3 次元のグラフを２次元にプロットするグラフです．上のコマン ドは，x，ｙの関数 f を x の値が xmin から xmax まで， ｙの値が ymin から ymax までの範囲で計算し z の値の大 きさを濃淡で示すことで２次元グラフをかきます．‘（ ）’内は作図のオプションで省略できることを示します． (例) z=sinx*cosx を-< x < ，-< y < の範囲でプロットす る． DensityPlot[Sin[x]*Cos[y],{x,-Pi,Pi},{y,-Pi,Pi}] -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 等高線グラフ ContourPlot［ｆ，｛ｘ，ｘｍｉｎ，ｘｍａｘ｝,｛ｙ，ymin，ｙｍａｘ｝（，ｏｐｔｉｏｎ -> ｖａｌｕｅ）］ 等高線グラフは密度グラフに似ていますが，値の大きさが同じ点を結んだ線を引くことで 3 次元のグラフを２次 元にプロットするグラフです．上のコマンドは，x，ｙの関数 f を x の値が xmin から xmax まで， y の値が ymin か ら ymax までの範囲で計算し z=f(x,y)の値の大きさが同じ点を結んだ線を引き２次元グラフをかきます．‘（ ）’ 内は作図のオプションで省略できることを示します． (例) z=sinx*cosx を-< x < ，-< y < の範囲でプロットす る． ContourPlot[Sin[x]*Cos[y],{x,-Pi,Pi},{y,-Pi,Pi}] -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3
（問題 5） (x2 y2)
e
の３次元グラフを-2 < x < 2，-2 < y < 2 の範囲で書きなさい．また等高線図を書きなさい． （応用例） 点電荷の作るポテンシャルと電場 以下では，原点に１C の電荷を置いたときにできる x-y 面内のポテンシャルと電場の図を作成する例を示す．簡 単のために，4
0
1
とする． ポ テンシャルは，r
z
y
x
1
4
1
)
,
,
(
0
で与えられるので， x-y 面 内 で は ， 2 21
)
,
,
(
y
x
z
y
x
と な る （ 但 し ，1
4
0
とした）．これを
8
x
8
，
8
y
8
で 3 次元 でプロットする. Phi1=1/Sqrt[x^2+y^2] Plot3D[Phi1,{x,-8,8},{y,-8,8}] このままでは，値域が狭いので0
2
でプロットすること にしよう． Plot3D[Phi1,{x,-8,8},{y,-8,8},PlotRange->{0,2}] さらに，メッシュが荒いので，x-y 面を 40x40 で切って表示す ることにする． Plot3D[Phi1,{x,-8,8},{y,-8,8},PlotRange->{0,2}, PlotPoints->40] 等高線図を書いてみよう． ContourPlot[Phi1,{x,-8,8},{y,-8,8}, PlotRange->{0,2},PlotPoints->40] -5 0 5 -5 0 5 0.1 0.2 0.3 0.4 -5 0 5 -5 0 5 -5 0 5 0 0.5 1 1.5 2 -5 0 5 -5 0 5 -5 0 5 0 0.5 1 1.5 2 -5 0 5 -7.5 -5 -2.5 0 2.5 5 7.5 -7.5 -5 -2.5 0 2.5 5 7.5等高線が見やすいように色塗りをやめる．後で電場と重ねて 描くために，「g1」という名前をつけておく． g1=ContourPlot[Phi1,{x,-8,8},{y,-8,8}, PlotRange->{0,2},PlotPoints->40, ContourShading->False] 電場を書いてみよう．電場 E は，
E
grad
で与えられ る．Mathematica では関数の勾配を描画するために， GradientFieldPlot[] という関数が用意されている．これを使うためには最初に一 度だけ， <<VectorFieldPlots` とする必要がある．ここで[ ` ]はバッククウォートでクウォート [ ‘ ]（７の上のキー）ではありません．バッククウォートは「シ フトキーを押しながら＠キーを押して」入力します． (注) 「General::obspkg: VectorFieldPlots`はサポートされなくなりました．ロ ードしようとしているレガシーバージョンは，現在の Mathematica 機能と衝突を起 す可能性があります．更新情報については Compatibility Guide をご覧くださ い． >>」というエラーメッセージが出ますが，新しいバージョンのバグなので気に しないで下さい． ここではメッシュを 20x20 に切って勾配をプロットすることにす る． g2=GradientFieldPlot[-Phi1,{x,-8,8},{y,-8,8}, ScaleFunction->(.6&),PlotPoints->20] 最後に，ポテンシャルと電場を重ねて描いてみよう． Show[g1,g2,PlotRange->{{-8,8},{-8,8}}] -7.5 -5 -2.5 0 2.5 5 7.5 -7.5 -5 -2.5 0 2.5 5 7.5 -6 -4 -2 0 2 4 6 8 -6 -4 -2 0 2 4 6 8 （問題６） 上の例に倣って，x
2
とx
2
にそれぞれ１C の電荷を置いたときの x-y 面内でのポテンシャル の等高線と電場を重ねて描きなさい．但し4
0
1
とする． （問題７） 上の例に倣って，x
2
とx
2
にそれぞれ-1C と１C の電荷を置いたときの x-y 面内でのポテンシ ャルの等高線と電場を重ねて描きなさい．但し4
0
1
とする．8.3.4. 式の変換
式の展開 Expand［expr］ 式 expr を展開する． (例) 2)
1
(
x
を展開する Expand[(x+1)^2] 1 + 2 x + x２ 因数分解 Factor[expr] 式 expr を因数分解する． (例) 22
1
x
x
を因数分解する Factor[1 + 2 x + x^2] ( 1 + x )2 式の簡約 Simplify[expr] 式 expr をなるべく簡単な形に変換する． (例) 22
1
x
x
を簡単にする Simplify[1 + 2 x + x^2] ( 1 + x )2 通分 Together[expr] 分数式の和 expr を通分する (例)1
/
x
1
/
y
を通分する Together［1/x + 1/y］ x y x y 部分分数展開 Apart[expr] 分数式 expr を部分分数に変換する (例)xy
y
x
y
x
1
2
3
を部分分数に展開する Apart[(3+2x+y)/(1+x+y+x y)] 1 1 x 2 1 y （問題 8）次の式を因数分解しないさい （1） 4 2 2 4b
b
a
a
(2)a
3
b
3
c
3
3
abc
(3) 3 3 3 3)
(
a
b
c
a
b
c
8.3.5. 関数の極限
Limit[f,x->x0] は x の関数 f の x->x0 の極限f
x xlim
0 を求めます． （例） x xxe
lim
を求める． Limit[x*Exp[-x],x->Infinity] 0 (問題 9)x
x
xsin
lim
0 を求めなさい． (問題 10) x xx
1
1
lim
を求めなさい．8.3.6. テーラー展開
Series[f(x),{x,x0,n}] は x の関数 f(x)の x=x0 の周りのテーラー展開を n 次まで求めます． (例) sinx を x=0 のまわりに 5 次まで展開し，その結果をグラフにして sinx のグラフと比較する． In[1]:= u=Series[Sin[x],{x,0,5}] Out[1]= x x3 6 x5 120 O x 6 In[2]:= u=Normal[u] Out[2]= x -x3 6 + x5 120 In[3]:=Plot[{Sin[x],u},{x,0,Pi}] Out[3]= -Graphics- （問題 11） xe
の x=0 のまわりでのテイラー展開を１次から３次まで求め， xe
と一緒にプロットし 展開次数が上がれば近似が良くなることを確認しなさい．8.3.7. 級数和
Sum[ｆ，{i,imin.imax}]級数 f(i)を i が imin から imax までの合計をもとめ る． (例) 自然数を１から１０まで合計する． Sum[i,{i,1,10}] 55 （問題 12）
11
i ni
を n=1,2,3,4 についてもとめ，n=1 の時だけ発散することを確認しなさい （問題 13）
0 n nx
を求めなさい8.3.8. 微分・積分
偏微分 D[f,x] 関数 f の x についての偏微分を求める． (例)sin
x
の x についての偏微分 D[Sin[x],x] Cos[x] 高階偏微分 D[f,{x,n}] 関数 f の x についての n 階微分 (例) nx
の x についての 3 階微分 D[x^n,{x,3}] 2 n 1 n n x 3 n 多重偏微分 D[f,x,y,...] 関数 f の x，y，．．．についての多重偏微分 (例)
sin
xy
の x，y についての多重偏微分 D[Sin[x*y],x,y] Cos[x y] - x y Sin[x y] 不定積分 Integrate[f,x] 関数 f の x についての不定積分 (例）sin
x
のxについての不定積分 Integrate[Sin[x],x] -Cos[x] 定積分 Integrate[f,{x,xmin,xmax}] 関数 f の xmin から xmax までの定積分 (例)sin
x
の0
x
での定積分 Integrate[Sin[x]^2,{x,0,Pi}]2
数値定積分 NIntegrate[f,{x,xmin,xmax}] 関数 f の xmin から xmax までの数値定積分 (例)sin
2x
の0
x
での定積分 NIntegrate[Sin[x]^2,{x,0,Pi}] 1.5708 （問題１４）関数 2 21
)
(
y
x
x
f
について次の偏微分を求めなさい． (1)f
x (2)f
y (3)f
xx (4)f
xy (5)f
yy （問題１５）以下の不定積分を求めなさい． (1)
log
xdx
(2)
e
xsin
xdx
(3)
x
x
2
1
dx
（問題１６）以下の定積分を求めなさい (1)
0 2dx
e
x (2)dx
x
x
sin
8.3.9. 方程式の解
1 元方程式 Solve[{lhs == rhs},x] は lhs=rhs という方程式を変数 x について解きま す． 「==」は「=」を２つ書くことに注意 (例)x
2
x
1
0
を x について解く． Solve[x^2+x+1==0,x] x 1 1 3 , x 1 2 3 多元連立方程式 Solve[{lhs1 == rhs1,lhs2==rhs2,...},{x,y,...}] は lhs1=rhs1，lhs2=rhs2，．．．という x，y，．．． についての連立方程式を解きます． (例) 連立方程式
3
5
3
x
y
x
y
を解く Solve[{3 x+y==5,x+y==3},{x,y}] {{x -> 1, y -> 2}} （問題１７）質量 m1の質点が速さ v1で，静止している質量 m2の質点に完全弾性衝突し，それぞれ v2と v3 の速さで動き始めた時，運動量保存則と運動エネルギー保存則から，
2 3 2 2 2 1 2 1 1 3 2 2 1 1 12
1
2
1
2
1
v
m
v
m
v
m
v
m
v
m
v
m
が成り立つ．この方程式を解いて v2と v3を求めなさい。8.3.10. 微分方程式
DSolve[eqns,y[x],x] は x の関数 y についての微分方程式 eqns を解きます． 単振動の例 (例１)微分方程式x
'
'
x
(
t
)
を解く． 「''」はダブルクオーテーションではなく，「'」 を２つ書く (例 2)微分方程式x
'
'
x
(
t
)
を初期条件x
(
0
)
1
,
1
)
0
(
'
x
のもとに解いて答えを u に代入する． 結果を0
x
2
でプロットする． DSolve[x''[t]==-x[t],x[t],t] {x[t] -> C[1] Cos[t] + C[2] Sin[t]} u=DSolve[{x''[t]==-x[t],x[0]==1,x'[0]==1}, ｘ[t],t] {{x[t] -> Cos[t] + Sin[t]}} Plot[Evaluate[x[t]]/.u, {t,0,2 Pi}] 1 2 3 4 5 6 -1 -0.5 0.5 1( 例 3) 連 立 微 分 方 程 式
y
'
(
t
)
x
(
t
)
，)
(
)
(
'
t
y
t
x
)を初期条件y
(
0
)
0
，x
(
0
)
1
のも とに解く DSolve[{y'[t]==x[t],x'[t]==-y[t],x[0]==1, y[0]==0},{x[t],y[t]},t]{{x[t] -> Cos[t], y[t] -> Sin[t]}}
（問題１８）
x
'
'
(
t
)
x
'
(
t
)
x
(
t
)
を初期条件x
(
0
)
1
，x
'
(
0
)
0
の下で解き，結果を0
t
4
でプ ロットしなさい．グラフが全て見えるように工夫すること．（問題１９）
