Mathematica入門

Mathematica 入 門
Ｍａｔｈｅｍａｔｉｃａ入門
０．目次
１．はじめに
２．数値計算
２．１
四則演算と定数
２．２
組み込み関数
２．３
数の表示
２．４
代入
２．５
Ｐｒｉｎｔ文
２．６
Ｄｏ文
２．７
Ｆｏｒ文
２．８
Ｗｈｉｌｅ文
２．９
Ｉｆ文
２．１０ 和と積
２．１１ 関数の定義
［ 例 ２ ． １ ］ 1+1/2+1/3+… +1/n+… の 計 算
［例２．２］πの計算
［例２．３］２＾ｎ−１の素数
３．ファイル処理
３．１ ファイルからコマンド入力
３．２ ファイルへ結果の出力
［例３．１］ファイル処理の例
４．数式処理
４．１
４．２
４．３
４．４
方程式
式の展開、因数分解、部分分数展開、簡約化
級数、微分、定積分
漸化式
５．リスト処理
５．１ リスト
５．２ 素数
５．３ 集合
［例５．１］素数に関する問題
［例５．２］ピタゴラス数
［例５．３］リスト処理
６．グラフィックス
６．１ 関数のグラフ
６．２ データのグラフ
６．３ データファイルの読込み
［例６．１］グラフィックスの例
７．課題
７．１
７．２
７．３
７．４
課題１
課題２
課題３
課題４
- 1 -
Mathematica 入 門
１．はじめに
Mathematicaは 、 数 式 計 算 、 数 値 計 算 、 グ ラ フ ィ ッ ク ス 処 理 等 を 実 行 す る 数 式 処
理支援ソフトウェアである。例えば、分数での計算、因数分解、方程式の解、微
積分等が誤差なく厳密に計算できる。もちろん、数値計算も行える。簡単な使い
方を示す。
% math （ 起 動 ）
In[1]:= 1+2*3
Out[1]= 7
In[2]:= 1/2+1/3
5
Out[2]= 6
In[3]:= Sqrt[8]
Out[3]= 2 Sqrt[2]
In[4]:= N[Sqrt[8]]
Out[4]= 2.82843
In[5]:= %4/2
Out[5]= 1.41421
In[6]:= Solve[x^2-2==0,x]
Out[6]= {{x -> -Sqrt[2]}, {x -> Sqrt[2]}}
In[7]:=Quit
%
（停止）
① In[k]は 、 起 動 時 か ら k番 目 に 入 力 し た 数 式 を 意 味 し 、 Out[k]は 、 In[k]に
対 応 し た 出 力 を 意 味 す る 。 Mathematicaを 起 動 し て か ら の 累 計 番 号 （ kは 1か ら
始まる）で自動的に表示される。
② In[1]:= と 表 示 さ れ て い る と き 、 1+2*3 と 入 力 し 、 Enterキ ー を 押 す 。
Out[1]= 7 と 表 示 さ れ る 。
③ Out[k]の 結 果 を 再 利 用 す る と き は 、 %k と す る 。 %は 直 前 の 結 果 を 意 味 す る 。
④ 計 算 の 中 断 は 、 C-c（ Ctrlキ ー を 押 し な が ら cを 押 す ） を 入 力 す る 。
Interrupt>と 表 示 さ れ た ら 、 aを 入 力 す る 。
⑤式の最後にセミコロンをつけると、計算はするが結果は出力されない。
- 2 -
Mathematica 入 門
２．数値計算
数値計算に必要な基本的書き方、実行文を示す。
２．１
四則演算と定数
加法
減法
乗法
除法
ベキ乗
階乗
ａ＋ｂ
ａ−ｂ
ａ＊ｂ
ａ／ｂ
ａ＾ｂ
ｎ！
Pi
E
（または、ａ
ｂ）
% math
In[1]:= 123 + 45
Out[1]= 168
In[2]:= 1/2 + 1/3
5
Out[2]= 6
In[3]:= 123 - 45
Out[3]= 78
In[4]:= 1/2 - 1/3
1
Out[4]= 6
In[5]:= 123 * 45
Out[5]= 5535
In[6]:= (1/2) * (1/3)
1
Out[6]= 6
In[7]:= 123 / 45
41
Out[7]= -15
In[8]:= (1/2) / (1/3)
3
Out[8]= 2
In[9]:= 2 ^ 10
Out[9]= 1024
In[10]:= 10 !
Out[10]= 3628800
In[11]:= Pi
Out[11]= Pi
In[12]:= E
Out[12]= E
In[13]:= Quit
- 3 -
円周率
自然対数の底
Mathematica 入 門
２．２
組み込み関数
絶対値
剰余
平方根
指数関数
三角関数
対数関数
２項係数
Abs[x]
Mod[x,k] 整 数 ｘ を 整 数 ｋ で 割 っ た 余 り
Sqrt[x]
Exp[x]
Sin[x] Cos[x] Tan[x] ArcTan[x]
Log[x]（ 底 は e） Log[a,x]（ 底 は a）
Binomial[n,r]
% math
In[1]:= Abs [123]
Out[1]= 123
In[2]:= Abs [-123]
Out[2]= 123
In[3]:= Mod [10,3]
Out[3]= 1
In[4]:= Sqrt [2]
Out[4]= Sqrt[2]
In[5]:= Exp [1]
Out[5]= E
In[6]:= Sin [Pi]
Out[6]= 0
In[7]:= Cos [Pi]
Out[7]= -1
In[8]:= Tan [Pi/4]
Out[8]= 1
In[9]:= Log [E]
Out[9]= 1
In[10]:= Log [2,4]
Out[10]= 2
In[11]:= Binomial [10,5]
Out[11]= 252
In[12]:= Quit
２．３
数の表示
N[x]
N[x,k]
x//N
ｘを小数表示にする。
ｘをｋ桁まで正しく小数表示する。
ｘを仮数部と指数部で表示する。
In[1]:=
Out[1]=
In[2]:=
Out[2]=
In[3]:=
N [Sqrt[2]]
1.41421
N [Sqrt[2],10]
1.414213562
2^100 //N
30
Out[3]= 1.26765 10
- 4 -
Mathematica 入 門
２．４
代入
x=値
Clear[x]
In[1]:=
Out[1]=
In[2]:=
Out[2]=
In[3]:=
In[4]:=
Out[4]=
２．５
x = 123
123
x
123
Clear [x]
x
x
Ｐｒｉｎｔ文
Print[S]
In[1]:=
Out[1]=
In[2]:=
a bc
In[3]:=
Out[3]=
In[4]:=
123
２．６
変 数 xに 値 を 代 入 す る 。 後 で 参 照 で き る 。
変 数 xに 割 り 当 て ら れ て い る 値 を 消 去 。
文字列Ｓが出力される。
空 白 は 、 二 重 引 用 符 （ "） で 必 要 な 空 白 文 字 を 囲 む 。
x = "a bc"
a bc
Print [x]
x = 123
123
Print [x]
Ｄｏ文
Do[X(n),{n,a,b,c}]
a≦ n≦ bに つ い て 、 増 分 cで X(n)を 評 価 し た
結果を出力する。
In[1]:= Do [Print[i],{i,1,3}]
1
2
3
２．７
Ｆｏｒ文
For[i=a,i<=b,i=i+c,X(i)]
a≦ i≦ bに つ い て 、 増 分 cで X(i)を 評 価
した結果を出力する。
In[1]:= For [i=1,i<=3,i=i+1,Print[i]]
1
2
3
- 5 -
Mathematica 入 門
２．８
Ｗｈｉｌｅ文
While[条 件 ,処 理 ]
条件が真の間、処理を繰り返す。
In[1]:=i = 1
Out[1]= 1
In[2]:= While [i<=3,{Print[i]; i=i+1}]
1
2
3
（ 注 意 ） 処 理 が 複 数 の 文 か ら 構 成 さ れ る と き 、 {と }で 囲 む 。
条件を関係演算子、論理演算子で表す。
関係演算子
x == y
xと yが 等 し い
x != y
xと yが 等 し く な い
x < y
xが yよ り 小 さ い
x <= y
xが yよ り 小 さ い か 等 し い
x > y
xが yよ り 大 き い
x >= y
xが yよ り 大 き い か 等 し い
２．９
論理演算子
x ¦¦ y
xま た は y
x && y
xか つ y
!x
xで な い
Ｉｆ文
If[条 件 ,処 理 １ ,処 理 ２ ]
条件が真ならば処理１を実行し、偽な
らば処理２を実行する。
In[1]:= If [1<2,Print["1"],Print["2"]]
1
In[2]:= If [1>2,Print["1"],Print["2"]]
2
２．１０
和と積
Sum[f[n],{n,a,b}]
Product[f[n],{n,a,b}]
In[1]:=
Out[1]=
In[2]:=
Out[2]=
f(a)+f(a+1)+・ ・ ・ +f(b-1)+f(b)
f(a)・ f(a+1)・ ・ ・ f(b-1)・ f(b)
Sum [n,{n,1,10}]
55
Product [n,{n,1,10}]
3628800
- 6 -
Mathematica 入 門
２．１１
関数の定義
関 数 名 [引 数 _]:=定 義 式
関 数 名 [引 数 １ _,引 数 ２ _]:=定 義 式 （ 引 数 が ２ 個 以 上 の 場 合 も 同 様 ）
?関 数 名
関数定義の確認。
Clear[関 数 名 ]
関数の定義削除。
関 数 名 [引 数 _]:=定 義 式 /; 条 件
条件が満たされたとき定義される
関 数 名 [引 数 _]:= (コ マ ン ド ; … ;
コマンド)
複数のコマンドをまとめてひとつ
の手続きとする。丸かっこで囲む
関数の値は、最後のコマンドの結
果となる。
関 数 名 [引 数 _]:= Module[{u,v},
コ マ ン ド ; … ;コ マ ン ド ]
手 続 き 内 で 使 う 局 所 変 数 u,vを 定 義
できる。
Return[値 ]
値を返す。
% math
In[1]:= f[x_,y_]:= x*x + y*y
In[2]:= ?f
Global`f
f[x_, y_] := x*x + y*y
In[3]:= f[3,4]
Out[3]= 25
In[4]:= Clear [f]
In[5]:= ?f
Global`f
In[6]:= f[x_]:= 1 /; x>=0
In[7]:= ?f
Global`f
f[x_] := 1 /; x >= 0
In[8]:= f[3]
Out[8]= 1
In[9]:= f[-3]
Out[9]= f[-3]
（続く）
- 7 -
Mathematica 入 門
In[10]:= d[x_,y_]:=(z=x*x+y*y; z=Sqrt[z])
In[11]:= d[2,3]
Out[11]= Sqrt[13]
In[12]:= a = 3
Out[12]= 3
In[13]:= e[x_]:=Module[{u},u=a*a+x*x; a=-a; u=Sqrt[u]]
In[14]:= e[4]
Out[14]= 5
In[15]:= a
Out[15]= -3
In[16]:= u
Out[16]= u
In[17]:= Quit
（参考）ヘルプ機能
？関数名
？？関数名
？文字列＊
？＊文字列＊
関数の簡単な説明
関数の詳しい説明
文字列で始まる関数の一覧
文字列を含む関数の一覧
- 8 -
Mathematica 入 門
［ 例 ２ ． １ ］ 1+1/2+1/3+… +1/n+… の 計 算
% math
Mathematica 3.0 for Digital Unix
Copyright 1988-97 Wolfram Research, Inc.
License valid through 24 May 2003.
-- Terminal graphics initialized -一般項の定義
In[1]:= a[k_]:=1/k;
In[2]:= b[n_]:=Sum[a[k],{k,1,n}];
In[3]:= a[3]
1
Out[3]= 3
In[4]:= b[3]
11
Out[4]= -6
分数表現
In[5]:= b[100]
14466636279520351160221518043104131447711
Out[5]= ----------------------------------------2788815009188499086581352357412492142272
小数表現
In[6]:= b[100]//N
Out[6]= 5.18738
In[7]:= b[1000]//N
Out[7]= 7.48547
In[8]:= b[10000]//N
Out[8]= 9.78761
In[9]:= Quit
- 9 -
Mathematica 入 門
［例２．２］πの計算
π /4 = 1 - 1/3 + 1/5 - 1/7 + 1/9 - 1/11 + 1/13 …
% math
Mathematica 3.0 for Digital Unix
Copyright 1988-97 Wolfram Research, Inc.
License valid through 8 Apr 1999.
-- Motif graphics initialized -一般項の定義。
In[1]:= f[k_]:=(-1)^(k+1)/(2*k-1);
５０番目の項までの部文和。
In[2]:= 4*Sum[f[k],{k,1,50}]
3400605476464206445954873476681150352328
Out[2]= ---------------------------------------1089380862964257455695840764614254743075
数値化。
In[3]:= N[%,20]
Out[3]= 3.1215946525910104785
In[4]:= g[n_]:=4*Sum[f[k],{k,1,n}];
１０番目のまでの項の部文和、２０番目までの項の部文和、・・・、５０番目ま
での項の部文和を求める。
In[5]:= Do[Print["g(",n,")=",N[g[n],20]],{n,10,50,10}]
g(10)=3.0418396189294022111
g(20)=3.0916238066678386317
g(30)=3.1082685666989461300
g(40)=3.1165965567938323178
g(50)=3.1215946525910104785
収束が遅いことがわかる。
In[6]:= Quit
- 10 -
Mathematica 入 門
［例２．３］２＾ｎ−１の素数
% math
Mathematica 3.0 for Digital Unix
Copyright 1988-97 Wolfram Research, Inc.
License valid through 24 May 2003.
-- Terminal graphics initialized -関 数 の 定 義 （ PrimeQ[x]は xが 素 数 な ら 真 、 素 数 で な い な ら 偽 と な る 関 数 ）
In[1]:= p[m_]:=For[n=1,n<=m,n++,
If[PrimeQ[2^n-1],Print["2^",n,"-1 is prime"]]];
関数の確認
In[2]:= ?p
Global`p
p[m_] := For[n = 1, n <= m, n++,
If[PrimeQ[2^n - 1], Print["2^", n,"-1 is prime"]]]
1000以 下 に つ い て 実 行 。
In[3]:= p[1000]
2^2-1 is prime
2^3-1 is prime
2^5-1 is prime
2^7-1 is prime
2^13-1 is prime
2^17-1 is prime
2^19-1 is prime
2^31-1 is prime
2^61-1 is prime
2^89-1 is prime
2^107-1 is prime
2^127-1 is prime
2^521-1 is prime
2^607-1 is prime
Timing[処 理 ]で 処 理 の 所 要 時 間 を 求 め る こ と が で き る 。
In[4]:= Timing[PrimeQ[2^607-1]]
Out[4]= {0.149994 Second, True}
In[5]:= Timing[PrimeQ[2^1000-1]]
Out[5]= {0. Second, False}
In[6]:= Quit
- 11 -
Mathematica 入 門
３．ファイル処理
ファイルからコマンドを入力する方法、ファイルへ実行結果を出力する方法を
示す。
３．１
ファイルからコマンド入力
フ ァ イ ル （ fib.m）
fib[0]:=1; (* コ メ ン ト *)
fib[1]:=1;
fib[n_]:=fib[n-1]+fib[n-2]/; n>1
こ の フ ァ イ ル を 読 み 込 む と 、 fibと い う 名 前 の 関 数 が 使 え る よ う に な る 。
(*と *)で 囲 ま れ た 部 分 は コ メ ン ト と 見 な さ れ 、 計 算 の 対 象 と な ら な い 。
!!fib.m
Get["fib.m"]
?fib
fib[10]
３．２
フ ァ イ ル （ fib.m） の 表 示 。
フ ァ イ ル （ fib.m） を 読 み 込 む 。
fibと い う 関 数 が 定 義 さ れ て い る か 確 認 す る 。
fib[10]を 計 算 す る 。
ファイルへ結果の出力
Put[式 ,"output"]
PutAppend[式 ,"output"]
Save["func.m",fib]
式 の 結 果 を 出 力 フ ァ イ ル (output） へ 保 存 。
式 の 結 果 を 出 力 フ ァ イ ル （ output） へ 追 加 保 存
フ ァ イ ル （ func.m） に 関 数 fibを 保 存 す る 。
こ の フ ァ イ ル を 読 み 込 め ば 、 関 数 fibが 使 え る 。
- 12 -
Mathematica 入 門
［例３．１］ファイル処理の例
% cat fib.m
fib[0]:=1;
fib[1]:=1;
fib[n_]:=fib[n-1]+fib[n-2]/; n>1
% math
フ ァ イ ル （ fib.m） を 読 み 込 む 。
In[1]:= Get["fib.m"]
関 数 fibが 定 義 さ れ て い る か 確 認 。
In[2]:= ?fib
Global`fib
fib[0] := 1
fib[1] := 1
fib[n_] := fib[n - 1] + fib[n - 2] /; n > 1
fib[5]を 計 算 。
In[3]:= fib[5]
Out[3]= 8
fib[5]の 結 果 を フ ァ イ ル （ output） に 出 力 。
In[4]:= Put[fib[5],"output"]
フ ァ イ ル （ output） を 表 示 。
In[5]:= !!output
8
fib[6]の 結 果 を フ ァ イ ル （ output） に 追 加 出 力 。
In[5]:= PutAppend[fib[6],"output"]
フ ァ イ ル （ output） を 表 示 。
In[6]:= !!output
8
13
フ ァ イ ル （ func.m） に 関 数 fibを 保 存 。
In[6]:= Save["func.m",fib]
フ ァ イ ル （ func.m） を 表 示 。
In[7]:= !!func.m
fib[0] := 1
fib[1] := 1
fib[n_] := fib[n - 1] + fib[n - 2] /; n > 1
In[7]:= Quit
- 13 -
Mathematica 入 門
４．数式処理
数式処理に必要な基本的な書き方、実行文を示す。
４．１
方程式
Solve[P==Q,x]
NSolve[P==Q,x]
N[Solve[P==Q,x],k]
Solve[{P1==Q1,P2==Q2},
{x1,x2}]
NSolve[{P1==Q1,P2==Q2},
{x1,x2}]
N[Solve[{P1==Q1,P2==Q2},
{x1,x2}],k]
FindRoot[P1==Q1,{x,a}]
FindMinimum[P,{x,a}]
P,Qを xに つ い て の 式 と す る 。
代数的な解が求められる。
P,Qを xに つ い て の 式 と す る 。
数値的な解が求められる。
P,Qを xに つ い て の 式 と す る 。
k桁 の 精 度 で 数 値 的 な 解 が 求 め ら れ る 。
P1,Q1,P2,Q2を x1,x2に つ い て の 式 と す る 。
代数的な解が求められる。
P1,Q1,P2,Q2を x1,x2に つ い て の 式 と す る 。
数値的な解が求められる。
P1,Q1,P2,Q2を x1,x2に つ い て の 式 と す る 。
k桁 の 精 度 で 数 値 的 な 解 が 求 め ら れ る 。
P1,Q1を xに つ い て の 式 と す る 。
x=aを 初 期 値 と し て 数 値 解 を 探 す 。
Pを xに つ い て の 式 と す る 。
x=aを 初 期 値 と し て 極 小 値 を 探 す 。
% math
In[1]:= Solve [ x^2 + x + 1 == 0, x]
1/3
2/3
Out[1]= {{x -> -(-1)
}, {x -> (-1)
}}
In[2]:= NSolve [ x^2 + x + 1 == 0, x]
Out[2]= {{x -> -0.5 - 0.866025 I}, {x -> -0.5 + 0.866025 I}}
In[3]:= Solve [{x^2 + y^2 == 4, x + y == 1},{x,y}]
1 - Sqrt[7]
1 + Sqrt[7]
Out[3]= {{x -> -----------, y -> -----------},
2
2
>
1 + Sqrt[7]
1 - Sqrt[7]
{x -> -----------, y -> -----------}}
2
2
In[4]:= NSolve [{x^2 + y^2 == 4, x + y == 1},{x,y}]
Out[4]= {{x -> -0.822876, y -> 1.82288},
{x -> 1.82288, y -> -0.822876}}
In[5]:= FindRoot [ Sin[x] == Cos[x],{ {x,0}]
Out[5]= {x -> 0.785398}
In[6]:= FindRoot [ Sin[x] == Cos[x], {x,Pi}]
Out[6]= {x -> 3.92699}
In[7]:= FindMinimum [x^3 - 5*x^2 - 2*x + 10, {x,0}]
Out[7]= {-15.3778, {x -> 3.52259}}
In[8]:= Quit
- 14 -
Mathematica 入 門
４．２
式の展開、因数分解、部分分数展開、簡約化
Expand[P]
Factor[P]
Apart[P]
Simplify[式 ]
多項式Ｐを展開する。
多項式Ｐを因数分解する。
多項式Ｐを部分分数に展開する。
式を最小の要素（演算記号、かっこ、文字）で表現
する。
% math
In[1]:= Expand [(x+y)^4]
4
Out[1]= x
3
2
+ 4 x
y + 6 x
2
y
3
+ 4 x y
4
+ y
In[2]:= Factor [x^8 - 1]
2
4
Out[2]= (-1 + x) (1 + x) (1 + x ) (1 + x )
In[3]:= Factor [a^3 + b^3 + c^3 - 3*a*b*c]
2
Out[3]= (a + b + c) (a
2
- a b + b
2
- a c - b c + c )
In[4]:= Apart [1/(1-x^2)]
-1
1
Out[4]= ---------- + --------2 (-1 + x)
2 (1 + x)
In[5]:= Simplify [(a+b)*(a-b)*(a-b)*(a+b)]
2
Out[5]= (a - b)
2
(a + b)
In[6]:= Simplify [Sin[x]*Cos[y]+Cos[x]*Sin[y]]
Out[6]= Sin[x + y]
In[7]:= Quit
- 15 -
Mathematica 入 門
４．３
級数、微分、微分方程式、定積分
Series[関 数 ,{x,a,n}]
関 数 を x=aの 近 傍 で n次 ま で べ き 級 数
展開する。
n階 導 関 数 を 求 め る 。
微分方程式を解く。
定積分を求める。
数値積分を求める。
D[関 数 ,{x,n}]
DSolve[微 分 方 程 式 ,関 数 ,変 数 ]
Integrate[関 数 ,{x,下 限 ,上 限 }]
NIntegrate[関 数 ,{x,下 限 ,上 限 }]
% math
In[1]:= Series [1/(1-2*x),{x,0,6} }]
2
Out[1]= 1 + 2 x + 4 x
3
+ 8 x
4
+ 16 x
5
6
7
+ 32 x
+ 64 x
+ O[x]
5
6
7
In[2]:= Series [1/(1-x-x^2),{x,0,6}]
2
Out[2]= 1 + x + 2 x
3
+ 3 x
4
+ 5 x
+ 8 x
In[3]:= Series [Sin[x],{x,0,6}]
3
5
x
x
7
Out[3]= x - -- + --- + O[x]
6
120
In[4]:= D [x^3,{x,1}]
2
Out[4]= 3 x
In[5]:= D [x^3,{x,2}]
Out[5]= 6 x
In[6]:= D [Log[x],{x,1}]
1
Out[6]= x
（続く）
- 16 -
+ 13 x
+ O[x]
Mathematica 入 門
In[7]:= DSolve [y'[x] == x^3,y[x],x]
4
x
Out[7]= {{y[x] -> -- + C[1]}}
4
In[8]:= DSolve [y'[x] == Cos[x],y[x],x]
Out[8]= {{y[x] -> C[1] + Sin[x]}}
In[9]:= Integrate [Sqrt[1-x^2],{x,0,1}]
Pi
Out[9]= -4
In[10]:= NIntegrate [Sqrt[1-x^2],{x,0,1}]
Out[10]= 0.785398
In[11]:= Quit
- 17 -
Mathematica 入 門
４．４
漸化式
Mathematicaに は 標 準 の 関 数 が 組 み 込 ま れ て い る が 、 パ ッ ケ ー ジ （ Mathematica
言語で書かれたプログラムファイルで、関数定義の集まりである）を読み込めば
機能を拡張できる。
RSolve[{初 期 値 ,漸 化 式 },関 数 ,変 数 ]
漸化式を初期値のもとで解く。
% math
In[1]:= Needs["DiscreteMath`RSolve`"]
In[2]:= RSolve[{f[1]==1,f[n]==f[n-1]+n},f[n],n]
n (1 + n)
Out[2]= {{f[n] -> ---------}}
2
In[3]:= RSolve[{f[1]==1,f[n]==2*f[n-1]+1},f[n],n]
n
Out[3]= {{f[n] -> If[n >= 1, -1 + 2 , 0]}}
In[4]:= RSolve[{f[0]==1,f[n]==n*f[n-1]},f[n],n]
Out[4]= {{f[n] -> n!}}
In[5]:= RSolve[{f[0]==1,f[1]==1,f[n]==f[n-1]+f[n-2]},f[n],n]
1 + n
1 - Sqrt[5] n
(1 + Sqrt[5])
(-----------) (-1 + Sqrt[5]) + -----------------2
n
2
Out[5]= {{f[n] -> --------------------------------------------------}}
2 Sqrt[5]
In[6]:= RSolve[{f[0]==1,f[1]==(1-Sqrt[5])/2,f[n]==f[n-1]+f[n-2]},f[n],n]
1
Sqrt[5] n
Out[6]= {{f[n] -> (- - -------) }}
2
2
In[7]:= Quit
- 18 -
Mathematica 入 門
５．リスト処理
５．１
リスト
数 、 変 数 、 式 な ど を コ ン マ （ ,） で 区 切 っ て 並 べ 、 ｛ ｝ で 囲 ん だ も の を リ ス ト
と い う 。 た と え ば 、 {2,3,5,7,11,13}, {x,y,z,1,2,3}, {{1,2},W,{a,b}}。
Part[L,i]
Length[L]
Sort[L]
リ ス ト Lの i番 目 の 要 素 。 L[[i]]と も 書 く 。
リ ス ト Lに 含 ま れ る 要 素 数 。
リストＬの要素を標準的順序に並べる。
% math
In[1]:= L = {11,22,33,44,55}
Out[1]= {11, 22, 33, 44, 55}
In[2]:= Part [L,3]
Out[2]= 33
In[3]:= L[[3]]
Out[3]= 33
In[4]:= Length [L]
Out[4]= 5
In[5]:= Sort [{3,5,2,1,4}]
Out[5]= {1, 2, 3, 4, 5}
In[6]:= Sort [{c,e,b,a,d}]
Out[6]= {a, b, c, d, e}
In[7]:= Sort [{a,23,2.3,aa}]
Out[7]= {2.3, 23, a, aa}
In[8]:= Sort [{{2,3},3,{1,2,3},{2,4}}]
Out[8]= {3, {2, 3}, {2, 4}, {1, 2, 3}}
- 19 -
Mathematica 入 門
Insert[L,x,n]
Delete[L,n]
First[L]
Take[L,n]
Rest[L]
Join[L1,L2,… ]
Union[L1,L2,… ]
リ ス ト Lの 前 か ら n番 目 に 要 素 xを 挿 入 す る 。
リ ス ト Lの 前 か ら n番 目 の 要 素 を 削 除 す る 。
リ ス ト Lの 先 頭 要 素 。
リ ス ト Lの 前 か ら n番 目 の 要 素 を 取 り 出 し リ ス ト に
する。
リ ス ト Lか ら 先 頭 要 素 を 除 い た 残 り の リ ス ト 。
リ ス ト L1,L2,… を こ の 順 に 結 合 す る 。
リ ス ト L1,L2,… を こ の 順 に 結 合 す る 。 た だ し 、 重
複する要素がある場合、１つにまとめ並べ替える。
In[9]:= M = {a,b,c,d,e}
Out[9]= {a, b, c, d, e}
In[10]:= Insert [M,f,1]
Out[10]= {f, a, b, c, d, e}
In[11]:= Delete [M,2]
Out[11]= {a, c, d, e}
In[12]:= First [{1,2,3}]
Out[12]= 1
In[13]:= Take [{1,2,3},1]
Out[13]= {1}
In[14]:= Rest [{1,2,3}]
Out[14]= {2, 3}
In[15]:= Rest [{{1,2},3,4}]
Out[15]= {3, 4}
In[16]:= Join [{2,3,4},{1,2,3}]
Out[16]= {2, 3, 4, 1, 2, 3}
In[17]:= Union [{2,3,4},{1,2,3}]
Out[17]= {1, 2, 3, 4}
In[18]:= Quit
- 20 -
Mathematica 入 門
Flatten[L]
Flatten[L,n]
Table[式 ,{i,a,b,d}]
Cases[L,形 式 ]
Cases[L,形 式 /;条 件 ]
リ ス ト Lの 入 れ 子 を 全 て 解 除 す る 。
リ ス ト Lの 入 れ 子 を 深 さ nま で 解 除 す る 。
iを aか ら bま で 増 分 dで 反 復 し 、 式 の 値 を リ ス ト に
する。
リ ス ト Lの 中 で 形 式 に 合 う も の を 選 択 す る 。
形式には、
x_
任 意 の も の （ 仮 に xと す る ）
x_Integer 整 数
x_List
リスト
などがある。
形 式 /;条 件 は 、 条 件 を 満 た す 形 式 を 意 味 す る 。
% math
In[1]:= Table [i^2, {i,1,9}]
Out[1]= {1, 4, 9, 16, 25, 36, 49, 64, 81}
In[2]:= Table [{i,i^2},{i,1,9}]
Out[2]= {{1, 1}, {2, 4}, {3, 9}, {4, 16}, {5, 25}, {6, 36}, {7, 49},
>
{8,64}, {9, 81}}
In[3]:= L = Table [{i,j},{1 i,1,2},{j,1,3}]
Out[3]= {{{1, 1}, {1, 2}, {1, 3}}, {{2, 1}, {2, 2}, {2, 3}}}
In[4]:= Flatten [L,1]
Out[4]= {{1, 1}, {1, 2}, {1, 3}, {2, 1}, {2, 2}, {2, 3}}
In[5]:= Flatten [L,2]
Out[5]= {1, 1, 1, 2, 1, 3, 2, 1, 2, 2, 2, 3}
In[6]:= Flatten [L]
Out[6]:= {1, 1, 1, 2, 1, 3, 2, 1, 2, 2, 2, 3}
In[7]:= M = Table [i,{i,1,10}]
Out[7]= {1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
In[8]:= Cases [M,x_/;Mod[x,2]==0]
Out[8]= {2, 4, 6, 8, 10}
In[9]:= Cases [{1,2,3.3,4,5.5},x_Integer]
Out[9]= {1, 2, 4}
In[10]:= Quit
- 21 -
Mathematica 入 門
５．２
素数
Prime[k]
PrimeQ[数 ]
GCD[i,j]
k番 目 の 素 数
数が素数か否か判定する。
素 数 な ら True、 合 成 数 な ら Falseを 返 す 。
ｉとｊの最大公約数
%math
In[1]:= Prime[1]
Out[1]= 2
In[2]:= Prime[10]
Out[2]= 29
In[3]:= PrimeQ[11]
Out[3]= True
In[4]:= PrimeQ[11111]
Out[4]= False
In[5]:= GCD[123,45]
Out[5]= 3
In[6]:= GCD[11,31]
Out[6]= 1
In[7]:= Quit
- 22 -
Mathematica 入 門
５．３
集合
集合はリストで表される。
Union[L1,L2]
Intersection[L1,L2]
Complement[S,L]
Outer[List,L1,L2]
Permutations[L]
リ ス ト L1と リ ス ト L2の 和 集 合 。
リ ス ト L1と リ ス ト L2の 共 通 部 分 。
全 体 集 合 Sに 関 す る 集 合 Lの 補 集 合 。
リ ス ト L1と リ ス ト L2の 直 積 。
リ ス ト Lの 要 素 の す べ て の 順 列
% math
In[1]:= Union[{1,2,3},{3,4,5}]
Out[1]= {1, 2, 3, 4, 5}
In[2]:= Intersection[{1,2,3,4},{3,4,5,6}]
Out[2]= {3, 4}
In[3]:= Complement[{1,2,3,4,5,6,7,8,9},{2,4,6,8}]
Out[3]= {1, 3, 5, 7, 9}
In[4]:= Outer[List,{a,b},{1,2}]
Out[4]= {{{a, 1}, {a, 2}}, {{b, 1}, {b, 2}}}
In[5]:= Permutations[{1,2,3}]
Out[6]= {{1, 2, 3}, {1, 3, 2}, {2, 1, 3}, {2, 3, 1}, {3, 1, 2},
{3, 2, 1}}
In[7]:= Quit
- 23 -
Mathematica 入 門
［例５．１］素数に関する問題
% math
１以上３０以下の素数を求める。
In[1]:= Cases[Table[i,{i,1,30}], x_/; PrimeQ[x]]
Out[1]= {2, 3, 5, 7, 11, 13, 17, 19, 23, 29}
1≦ x≦ 10で 、 x 2 +x+41 が 素 数 と な る xを 求 め る 。
In[2]:= Cases[Table[i,{i,1,10}], x_/; PrimeQ[x^2+x+41]]
Out[2]= {1, 2, 3, 4, 5, 6, 7, 8, 9, 10}
1≦ x≦ 100で 、 x 2 +x+41が 素 数 と な ら な い xを 求 め る 。
In[3]:= Cases[Table[i,{i,1,100}], x_/; !PrimeQ[x^2+x+41]]
Out[3]= {40, 41, 44, 49, 56, 65, 76, 81, 82, 84, 87, 89, 91, 96}
1≦ x≦ 100で 、 x 2 -79x+1601が 素 数 と な ら な い xを 求 め る 。
In[4]:= Cases[Table[i,{i,1,100}], x_/; !PrimeQ[x^2-79*x+1601]]
Out[4]= {80, 81, 84, 89, 96}
In[5]:= Quit
%
- 24 -
Mathematica 入 門
［例５．２］ピタゴラス数
a,b,cを 正 整 数 と す る 。 a,b,cが 、 a 2 +b 2 =c 2 (a≦ b≦ c)を 満 た す と き 、 ピ タ ゴ ラ ス
数 と い う 。 正 整 数 nが 与 え ら れ た と き 、 c 2 ≦ nを 満 た す ピ タ ゴ ラ ス 数 を 求 め よ 。
●方法１
% math
n=900の 場 合 。
３ つ 組 (a,b,c)の 集 合 L1を 求 め る 。 た だ し 、 1≦ a≦ 30,1≦ b≦ 30,1≦ c≦ 30。
In[1]:= L1=Flatten[Table[{i,j,k},{i,1,30},{j,1,30},{k,1,30}],2];
集 合 L1の 中 で 、 条 件 ： a 2 +b 2 =c 2 を 満 た す ３ つ 組 を 求 め る 。
In[2]:= Cases[L1, x_/; x[[1]]^2+x[[2]]^2==x[[3]]^2]
Out[2]= {{3, 4, 5}, {4, 3, 5}, {5, 12, 13}, {6, 8, 10}, {7, 24, 25},
>
{8, 6, 10}, {8, 15, 17}, {9, 12, 15}, {10, 24, 26}, {12, 5, 13},
>
{12, 9, 15}, {12, 16, 20}, {15, 8, 17}, {15, 20, 25}, {16, 12, 20},
>
{18, 24, 30}, {20, 15, 25}, {20, 21, 29}, {21, 20, 29}, {24, 7, 25},
>
{24, 10, 26}, {24, 18, 30}}
集 合 L1の 中 で 、 条 件 ： a 2 +b 2 =c 2 を 満 た す ３ つ 組 の 集 合 L2を 求 め る 。
In[3]:= L2= Cases[L1, x_/; x[[1]]^2+x[[2]]^2==x[[3]]^2];
集 合 L2の 中 で 、 条 件 ： (a≦ b≦ c)を 満 た す ３ つ 組 の 集 合 を 求 め る 。
In[4]:= Cases[L2, x_/; (x[[1]]<=x[[2]])&&(x[[2]]<=x[[3]])]
Out[4]= {{3, 4, 5}, {5, 12, 13}, {6, 8, 10}, {7, 24, 25}, {8, 15, 17},
>
{9, 12, 15}, {10, 24, 26}, {12, 16, 20}, {15, 20, 25}, {18, 24, 30},
>
{20, 21, 29}}
In[5]:= Quit
- 25 -
Mathematica 入 門
●方法２
(* pita1.m *)
pita1[n_] := (
count = 0;
For[c=1, c*c<=n, c++, {
For[a=1, 2*a*a<=c*c, a++ {
For[b=a, b<=c, b++, {
If[a^2+b^2==c^2,{ count++; }]
}];
}];
}];
Print["n=",n," count=",count];
)
実行結果
% math
In[1]:= Get["pita1.m] "]
In[2]:= pita1[100]
n=100 count=2
In[3]:= pita1[1000]
n=1000 count=11
In[4]:= pita1[10000]
n=10000 count=52
In[5]:= Timing[pita1[10000]]
n=10000 count=52
Out[5]= {15.9327 Second, Null}
In[6]:= Quit
- 26 -
Mathematica 入 門
［例５．３］リスト処理
与えられたリストを反転する関数を定義せよ。
{a, b, c, d} ---> {d, c, b, a}
●方法１
(* list11.m *)
f[L_] := (
If[ L == {},
{},
Join[f[Rest[L]],Take[L,1]]
]
)
実行結果
% math
In[1]:= Get["list11.m"]
In[2]:= f[{a}]
Out[2]= {a}
In[3]:= f[{a,b}]
Out[3]= {b, a}
In[4]:= f[{a,b,c}]
Out[4]= {c, b, a}
In[5]:= f[{a,b,c,d}]
Out[5]= {d, c, b, a}
In[6]:= Quit
- 27 -
Mathematica 入 門
●方法２
(* list12.m *)
g[L_] := Module[{W},
If[ L == {},
{},
{ W = {};
For[i=1, i<=Length[L], i++, {
W = Insert[W,Part[L,i],1];
}];
}
];
Return[W];
]
実行結果
% math
In[1]:= Get["list12.m"]
In[2]:= g[{a}]
Out[2]= {a}
In[3]:= g[{a,b}]
Out[3]= {b, a}
In[4]:= g[{a,b,c}]
Out[4]= {c, b, a}
In[5]:= Quit
- 28 -
Mathematica 入 門
６．グラフィックス
６．１
関数のグラフ
Plot[f[x],{x,a,b}]
Plot[{f[x],g[x]},{x,a,b}]
F=Plot[f[x],{x,a,b}];
Show[F]
ParametricPlot[{x[t],y[t]},
{t,t0,t1}]
関 数 f[x]を 区 間 a≦ x≦ bで 描 く 。
関 数 f[x],g[x]を 区 間 a≦ x≦ bで 描 く 。
区 間 a≦ x≦ bの 関 数 f[x]を Fに 保 存 。
Fを 表 示 。
パ ラ メ ー タ 表 示 。 t0≦ t≦ t1。
（ 注 意 ） グ ラ フ を 新 た な ウ ィ ン ド ウ 内 に 表 示 す る に は 、 mathematicaが Motif
graphicsと し て 起 動 さ れ て い る 必 要 が あ る 。
hcsを リ モ ー ト ホ ス ト と し て 使 う と き は 、
hcs% setenv DISPLAY ロ ー カ ル ホ ス ト :0
を実行しておくこと。
% math
Mathematica 3.0 for Digital Unix
Copyright 1988-97 Wolfram Research, Inc.
License valid through 24 May 2003.
-- Motif graphics initialized --
６．２
データのグラフ
ListPlot[リ ス ト ]
ListPlot[リ ス ト ,
PlotJoined->True]
６．３
リストをプロットする。
リストをプロットし折れ線で結ぶ。
データファイルの読込み
ReadList["フ ァ イ ル 名 ",
ファイルからデータを読込み、リスト
Number]
に変換する。
ReadList["フ ァ イ ル 名 ",
ファイルからデータを読込み、行ごと
Number,
にリストに変換する。
RecordLists->True]
- 29 -
Mathematica 入 門
［例６．１］グラフィックスの例
% math
Mathematica 3.0 for Digital Unix
Copyright 1988-97 Wolfram Research, Inc.
License valid through 24 May 2003.
-- Motif graphics initialized -sin関 数 の プ ロ ッ ト
In[1]:= Plot[Sin[x],{x,0,2Pi}]
Out[1]= -Graphics円のプロット（パラメトリック）
In[2]:= ParametricPlot[{Cos[x],Sin[x]},{x,0,2Pi}]
Out[2]= -GraphicsIn[3]:= ListPlot[{1.0,1.2,1.6,2.2,3.0,4.0,5.2,6.6}]
Out[3]= -GraphicsIn[4]:= ListPlot[{1.0,1.2,1.6,2.2,3.0,4.0,5.2,6.6},PlotJoined->True]
Out[4]= -GraphicsIn[5]:= !!data
1.0 1.2
1.6 2.2
3.0 4.0
5.2 6.6
ファイルから読込み
In[5]:= ReadList["data",Number]
Out[5]= {1., 1.2, 1.6, 2.2, 3., 4., 5.2, 6.6}
In[6]:= ReadList["data",Number,RecordLists->True]
Out[6]= {{1., 1.2}, {1.6, 2.2}, {3., 4.}, {5.2, 6.6}}
In[7]:= Quit
- 30 -
Mathematica 入 門
７．課題
７．１
課題１
フィボナッチ数について考察する。
f(n) = f(n-1) + f(n-2) (n≧ 2)
f(0) = 1, f(1) = 1
a = (1 + sqrt(5))/2, b = (1 - sqrt(5))/2 と す る 。
f(0) = 1, f(1) = 1 の 場 合 の 一 般 項
f(n) = (a^(n+1) - b^(n+1))/sqrt(5)
（ １ ） f(0) = 1, f(1) = b の 場 合 の 一 般 項
f(n) = b^n
（ １ ） の 場 合 、 bの 有 効 桁 数 を ６ 桁 程 度 と す る と 、 数 値 計 算 の 結 果 と 理 論 結 果 と が
異なることを説明せよ。
=====厳 密 解 =====
=====数 値 解 =====
f(10)=0.0081306187557833487477241098899
f(10)=0.00813062
f(20)=0.000066106961351895970065516458134
f(20)=0.000066107
-7
-7
f(30)=5.3749049985557033441446053199 10
f(30)=5.37445 10
-9
-9
f(40)=4.3701303391810674621873016337 10
f(40)=-1.18878 10
-11
-7
f(50)=3.5531863700963434577568956788 10
f(50)=-6.83665 10
-13
f(60)=2.88896037434990853680051125445 10
f(60)=-0.0000840896
-15
f(70)=2.34890354044042507733164210695 10
f(70)=-0.0103423
-17
f(80)=1.90980391814308317408537097662 10
f(80)=-1.27202
-19
f(90)=1.55278875567226793364194299097 10
f(90)=-156.449
-21
f(100)=1.26251333806384294218902270604 10
f(100)=-19241.9
- 31 -
Mathematica 入 門
(* m111.m *)
f[n_] := (
Print["=====厳 密 解 ====="];
f0=1;
f1=(1-Sqrt[5])/2;
For[i=2,i<=n,i++,{
f2=f1+f0;
If[Mod[i,10]==0, Print["f(",i,")=",N[f2,30]]];
f0 = f1; f1 = f2;
}];
Print["=====数 値 解 ====="];
f0=1;
f1=(1-Sqrt[5])/2;
f1 = N[f1];
For[i=2,i<=n,i++,{
f2=f1+f0;
If[Mod[i,10]==0, Print["f(",i,")=",f2]];
f0 = f1; f1 = f2;
}];
)
- 32 -
Mathematica 入 門
●誤差発生検証
フ ァ イ ル （ m112.m）
n:=100;
b1:=(1-Sqrt[5])/2 - 0.1^15;
Print["b1=(1-Sqrt[5])/2 - 0.1^15"];
h0:=1;
h1:=b1;
Do[h2=h1+h0;h0=h1; h1=h2,{i,2,n}]
Print["h(",n,")=",h2]
Print[" "]
b2:=(1-Sqrt[5])/2 - 0.1^16;
Print["b2=(1-Sqrt[5])/2 - 0.1^16"];
h0:=1;
h1:=b2;
Do[h2=h1+h0;h0=h1; h1=h2,{i,2,n}]
Print["h(",n,")=",h2]
Print[" "]
b3:=(1-Sqrt[5])/2 - 0.1^17;
Print["b3=(1-Sqrt[5])/2 - 0.1^17"];
h0:=1;
h1:=b3;
Do[h2=h1+h0;h0=h1; h1=h2,{i,2,n}]
Print["h(",n,")=",h2]
Print[" "]
b4:=(1-Sqrt[5])/2 - 0.1^18;
Print["b4=(1-Sqrt[5])/2 - 0.1^18"];
h0:=1;
h1:=b4;
Do[h2=h1+h0;h0=h1; h1=h2,{i,2,n}]
Print["h(",n,")=",h2]
Print[" "]
b5:=(1-Sqrt[5])/2 - 0.1^19;
Print["b5=(1-Sqrt[5])/2 - 0.1^19"];
h0:=1;
h1:=b5;
Do[h2=h1+h0;h0=h1; h1=h2,{i,2,n}]
Print["h(",n,")=",h2]
Print[" "]
- 33 -
Mathematica 入 門
% math
In[1]:= Get["m112.m"]
b1=(1-Sqrt[5])/2 - 0.1^15
h(100)=-373184.
b2=(1-Sqrt[5])/2 - 0.1^16
h(100)=-58568.8
b3=(1-Sqrt[5])/2 - 0.1^17
h(100)=-19241.9
b4=(1-Sqrt[5])/2 - 0.1^18
h(100)=-19241.9
b5=(1-Sqrt[5])/2 - 0.1^19
h(100)=-19241.9
In[2]:= Quit
- 34 -
Mathematica 入 門
７．２
課題２
ピ タ ゴ ラ ス 数 ： a 2 +b 2 =c 2 (a≦ b≦ c)
ｎ
10 2
10 3
10 4
10 5
10 6
で 、 c 2 ≦ nを 満 た す 個 数 を 効 率 よ く 求 め よ 。
個数
2
11
52
220
881
(* pita2.m *)
pita2[n_] := (
c = 1;
count = 0;
For[c=1, c*c<=n, c++, {
a = 1;
b = c;
While[ a<=b, {
w = a*a + b*b - c*c;
If[ w>0, b-- ];
If[ w == 0, count++ ];
a++;
}];
}];
Print["n=",n," count=",count];
)
% math
In[1]:= Get["pita2.m"]
In[2]:= pita2[10000]
n=10000 count=52
In[3]:= pita2[100000]
n=100000 count=220
In[4]:= pita2[1000000]
n=1000000 count=881
In[5]:= Timing[pita2[10000]]
n=10000 count=52
Out[5]= {0.366652 Second, Null}
In[6]:= Quit
- 35 -
Mathematica 入 門
７．３
課題３
（ １ ） m,nが 与 え ら れ た と き 、 1/a+1/b=m/n (a≦ b)を 満 た す a,bの 組 を
すべて求めよ。
実行結果
1/3 + 1/6
1/4 + 1/4
m=1
n=2
c=2
(* << egypt21.m >> *)
(* 1/a + 1/b = m/n *)
m = 1; n = 2;
count = 0;
a = Floor[n/m] + 1; a1 = Floor[2*n/m];
While[ a <= a1, {
m1 = m*a-n; n1 = a*n;
If[ Mod[n1,m1] == 0, {
count++;
Print["1/",a," + 1/",n1/m1];
}];
a++;
}];
Print["m=",m," n=",n," c=",count];
- 36 -
Mathematica 入 門
（ ２ ） 2≦ n≦ 8,1≦ m≦ n-1の 組 (m,n)に つ い て 、 1/a+1/b=m/n (a≦ b)を 満 た す a,bの
個数を求めよ。
実行結果
m=1
m=1
m=2
m=1
m=2
m=3
m=1
m=2
m=3
m=4
m=1
m=2
m=3
m=4
m=5
n=2
n=3
n=3
n=4
n=4
n=4
n=5
n=5
n=5
n=5
n=6
n=6
n=6
n=6
n=6
count=2
count=2
count=2
count=3
count=2
count=1
count=2
count=2
count=1
count=0
count=5
count=2
count=2
count=2
count=1
m=1
m=2
m=3
m=4
m=5
m=6
m=1
m=2
m=3
m=4
m=5
m=6
m=7
n=7
n=7
n=7
n=7
n=7
n=7
n=8
n=8
n=8
n=8
n=8
n=8
n=8
count=2
count=2
count=0
count=1
count=0
count=0
count=4
count=3
count=2
count=2
count=1
count=1
count=0
(* << egypt22.m >> *)
(* 1/a + 1/b = m/n *)
n= 2; nmax = 8;
While[ n <= nmax, {
m = 1;
While[ m < n, {
count = 0;
a = Floor[n/m] + 1; a1 = Floor[2*n/m];
While[ a <= a1, {
m1 = m*a-n; n1 = a*n;
If[ Mod[n1,m1] == 0, { count++; }];
a++;
}];
Print["m=",m," n=",n," count=",count];
m++;
}];
n++;
}];
- 37 -
Mathematica 入 門
（ ３ ） 2≦ n≦ 20,1≦ m≦ n-1の 分 数 m/nに つ い て 、 ３ 個 以 下 の 単 位 分 数 で 表 せ な い
ものを求めよ。
実行結果
8/11,9/11,10/11,
12/13,
13/14,
15/16,
8/17,14/17,15/17,16/17,
9/19,14/19,15/19,17/19,18/19
（ ４ ） 2≦ n≦ 30,1≦ m≦ n-1の 分 数 m/nに つ い て 、 ４ 個 以 下 の 単 位 分 数 で 表 せ な い
ものを求めよ。
実行結果
16/17,
21/23,22/23,
27/29,28/29
（ ５ ） 2≦ n≦ 100,1≦ m≦ n-1の 分 数 m/nに つ い て 、 ５ 個 以 下 の 単 位 分 数 で 表 せ な い
ものを求めよ。
実行結果
77/79
- 38 -
Mathematica 入 門
７．４
課題４
次の曲線を表示せよ。
（１）アステロイド
定義
x = a*cos(t)**3, y = a*sin(t)**3
（２）リサージュ図形
定義
x = cos(a*t), y = sin(b*t)
（３）平面曲線１
定義
x = (a-b)*cos(t) + b*cos((a-b)*t/b)
y = (a-b)*sin(t) + b*sin((a-b)*t/b)
（４）平面曲線２
定義
x = (a+b)*cos(t) - b*cos((a+b)*t/b)
y = (a+b)*sin(t) - b*sin((a+b)*t/b)
（５）平面曲線３
定義
r = a*sin(t/n) n:正 整 数
x = a*sin(t/n)*cos(t)
y = a*sin(t/n)*sin(t)
（６）円の伸開線
定義
x = a*cos(t) + a*t*sin(t)
y = a*sin(t) - a*t*cos(t)
円の伸開線は、円柱に巻いた糸をたるませないでほどいて
いくとき糸の先端が描く曲線。
（７）正葉形
定義
r = a*sin(n*t)
nが 偶 数 の と き 2n個 の 葉 を 持 ち 、 nが 奇 数 の と き n個 の 葉 を 持 つ
x = r*cos(t) = a*sin(n*t)*cos(t)
y = r*sin(t) = a*sin(n*t)*sin(t)
（８）リマソンのか牛線
定義
r = a*cos(t) + b
x = (a*cos(t)+b)*cos(t)
y = (a*cos(t)+b)*sin(t)
- 39 -
Mathematica 入 門
（９）コンコイド
定義
r = a/cos(t) + b
x = (a/cos(t)+b)*cos(t)
y = (a/cos(t)+b)*sin(t)
（１０）対数らせん
定義
r = a**t (a>1)
x = (a**t)*cos(t)
y = (a**t)*sin(t)
（１１）アルキメデスのらせん
定義
r = a*t (a>0)
x = a*t*cos(t)
y = a*t*sin(t)
（１２）双曲らせん
定義
r*t = a (a>0)
x = (a/t)*cos(t)
y = (a/t)*sin(t)
（１３）リチュウス
定義
(r**2)*t = a (a>0)
x = sqrt(a/t)*cos(t)
y = sqrt(a/t)*sin(t)
（１４）デカルトの葉形
定義
x = 3*a*t/(1+t**3)
y = 3*a*t**2/(1+t**3)
(a>0)
（１５）サイクロイド
定義
x = a*(t-sin(t))
y = a*(1-cos(t))
x軸 上 を 転 が る 半 径 aの 円 周 上 に 固 定 さ れ た 点 の 描 く 曲 線 。
（１６）簡単なフーリエ級数
定義
y = (2/pi)(sin(x)+sin(3*x)/3+sin(5*x)/5+・ ・ ・ )
- 40 -