としては計算ができない

(1)

2016年度・後期・計算機数学基礎・第８回 1

● 講義資料

▼ 前回の講義内容のコードに関して

例えば

ℓ_2n = 2n

√ 1−√

1−(ℓn/n)²

2 , ℓ₆ = 3

によって定義された数列の値を計算するには,

l[n_] := n Sqrt[(1 - Sqrt[1 - (l[n/2]/(n/2))^2])/2];

l[6] = 3;

によって ℓ_n を定義して, Table[l[6*2^n], n, 0, 10]

として値を得る. これを,

l[2n_] := 2n Sqrt[(1 - Sqrt[1 - (l[n]/n)^2])/2];

としては計算ができない.

Mathematica で“x_”なる記述は,「パターン」と呼ばれ,そのパターンにマッチする

ときに,対応する文字の置き換えが行われる. 例えば, f[a_+1] := f[2 a]; f[x+1]

とすることはできる（ただし,x に値が入っていない場合）. しかし,x に具体的な値が割り当てられると, f[x+1] を評価する時に,先に x+1 が評価されるので,a_+1 にはマッチしなくなる.

▼ モンテカルロ法による円周率の計算

• [0,1] 上の100 個の乱数をリストとして得る：

RandomReal[1, 100]

• [0,1]² 上の100 個の乱数をリストとして生成する：

RandomReal[1, {100, 2}]

• [0,1]² 上の 5000 個の乱数をリストとして得て, その点の中で単位円の内部に入る個数を得る：

n = 5000; pairs = RandomReal[1, {n, 2}];

normed = Map[Norm, pairs];

number = Count[normed, _?(#<=1&)]

• それらを図で表す：

list = ListPlot[pairs, AspectRatio -> 1];

region = RegionPlot[x^2 + y^2 < 1, {x, 0, 1}, {y, 0, 1}];

Show[{region, list}]

Dec. 02, 2016, Version: 1.1 [email protected]

(2)

• リスト sample に入ったデータの平均と不偏分散を得る：

Mean[sample]

Variance[sample]

��

��^��^�

��

(3)

▼ 円周率の１６進小数展開の任意の桁を得る

π =

∑∞ i=0

1 16ⁱ

( 4

8i+ 1 − 2

8i+ 4 − 1

8i+ 5 − 1 8i+ 6

)

(1)

• π の１６進小数展開の d桁目から 20 桁程度を表示する：

hexPi[d_] := Module[{p, l}, p = Pi*16^(d - 1);

l = IntegerLength[IntegerPart[p]];

BaseForm[FractionalPart[N[p, l + 22]], 16]];

d Mathematica 計算近似計算

1 243f6a8885a308d313 243f6a86f7efee1066 50 82efa98ec4e6c89453 82efa98ebf91e9316b 100 c29b7c97c50dd3f84d c29b7c97c3948885ad 1000 349f1c09b075372c98 349f1c09b0711103b4 10000 68ac8fcfb8016cbdb9 68ac8fcfb8016203ee 100000 535ea16c406363a30c 535ea16c4063638790

● 実習内容

1. [0,1]² に値を持つ n 個の乱数を使って, モンテカルロ法で π の値を計算するため

の関数

singleApproxPi[n_]

を書きなさい. 次に, それを m 回の試行を行った結果をリストとして得るための関数

approxPi[m_, n_]

を書きなさい. これを利用して, (n, m) = (100,100), (1000,1000) での π が含まれる 95% 信頼区間を求めなさい.

2. 式 (1) を用いて, π の１６進小数展開の d 桁め以下 20桁程度を表示する関数 approxHexPi[d_]

を書きなさい. さらに, それを用いて, １６進展開の d 桁め, d= 100, 1000, 10000 を求めなさい.

3. n= 4, 5 についてDⁿ ={x∈Rⁿ :|x|<1}の体積をモンテカルロ法を用いて求めなさい.

4. 積分

∫ π/2 0

√sin(x) + 2 cos(x)dx をモンテカルロ法を用いて求めなさい.