● 講義資料
▼ 講義予定
• 自然対数の底の近似値の計算と浮動小数点演算の誤差（続き）
• 円周率の近似値の計算と浮動小数点演算の誤差
● 前回の講義のまとめ
▼ ニュートン法について（続き）
• ニュートン法は以下の「逐次近似」の一例と考えられる.
• F:R−→Rを C2-級関数とする. この時,
xn+1 =F(xn)
によって数列{xn}を定義して,写像F の固定点(α=F(α)を満たす点)を近似する方法を 逐次近似と呼ぶ.
• 関数f(x)に対するニュートン法xn+1=xn−f(xn)/f′(xn)は, F(x) =x−f(x)/f′(x)と 定義した逐次近似である.
• Rのある近傍U と正定数L <1が存在して,任意のx,y∈U に対して
|F(x)−F(y)| ≤L|x−y|
を満たすとき,F はU 上で縮小写像であると言う.
• C2-級関数F:R−→Rが,U ⊂R上で縮小写像であるとき,F は U にただひとつの固定点 を持つ（縮小写像の原理）.
証明は,x0∈U,xn+1=F(xn)で定義した点列{xn}がコーシー列となることを示し,その 収束先が固定点となることを示せばよい.
• C2-級関数F: R−→Rに対する逐次近似xn+1 =F(xn)の固定点αがあらかじめわかって いる時を考える. この時,|F′(α)|<1が成り立てば,αの近傍でF は縮小写像となる. 特に, αの十分近くに初期値x0 をとった逐次近似xn+1=F(xn)は,固定点αに収束する.
証明は,xの近傍で平均値の定理を用いればよい.
• 逐次近似xn+1=F(xn)がαに収束すると仮定する. この時,ǫn=|xn−α|は ǫn+1≤ |F′(α)|ǫn+1
2|F′′(α)|ǫ2n+O(ǫ3n) をみたす. 証明は,αでF を３次までテイラー展開して計算すればよい.
• ニュートン法の収束の速さは,f(x) = 0の解が単根のとき, ǫn+1≤Cǫ2n
をみたし,m 重根の時
ǫn+1≤ m−1 m ǫn
をみたす.
証明は,F(x) =x−f(x)/f′(x)と考え,上の結果を適用すればよい.
• この結果から,ニュートン法は,解が単根の時と重根の時とでは,その収束の様子が質的に異 ることがわかる.
• さらに,ニュートン法の収束の判定には,次の結果を用いることができる.
1. ある値αの近似列{xn}が|α−xn+1| ≤C|α−xn|2を満たしていると仮定する. この時, 任意の0< C <1 に対して,|xn+1−xn|< Cǫ|xn|が成り立つならば,|α−xn|< ǫ|α|
が成り立つ.
2. ある値 α の近似列 {xn} が, ある 0 < L < 1 が存在して, |α−xn+1| ≤ L|α−xn| を満たしていると仮定する. この時, |xn+1−xn| < (1−L)ǫ|xn| が成り立つならば,
|α−xn|< ǫ|α|が成り立つ.
• たとえば,ニュートン法の単根の場合には,収束の判定条件として,「|xn+1−xn|> ǫ|xn|が成 り立つ間,繰り返しを続ける」と書けば,そのループから脱出した時点では,|xn+1−xn|< ǫ|xn| が成り立っているので,近似値xn+1 は解αを相対誤差ǫで近似していると考えてよいこと がわかる.
▼ 自然対数の底
e
の近似値の計算について• 今回考えるのは,
an=
1 + 1 n
n
の値を計算することによって,自然対数の底eの近似値を求める.
• 結論から言えば,この方法では必要な精度（倍精度浮動小数点数でeを十分に近似している と考えられる精度）で近似値を得ることができない.
• 数学的には
an =e− e
2n+O(n−2) が成り立つ.
証明は, log(1 +x) のテイラー展開にx= 1/nを代入し, さらにその結果をexp(x) のテイ ラー展開に適用すればよい.
• 実際に計算できるan の値（それをan と書こう）は,ある（小さな）δ >0 が存在して,
|an−e| ∼e(1 + 1 2n+nδ) を満たす.
● 講義資料
▼ 自然対数の底の近似値を求める
以下の図は, (1 + 1/n)nを計算して,eの近似値を求めた例である. （eと近似値の値との絶対誤 差を表示している）
1e-08 1e-07 1e-06 1e-05 0.0001 0.001 0.01 0.1 1
1 100 10000 1e+06 1e+08 1e+10 1e+12 1e+14 1e+16
高速乗算・pow関数を利用して計算した結果
1e-10 1e-08 1e-06 0.0001 0.01 1 100
1 100 10000 1e+06 1e+08 1e+10 1e+12 1e+14 1e+16
(1+1/n)^n Floating Point Arithmetics Error by Rounding Mode nearest upward downward toward zero
丸めモードを変更して計算した結果
▼ 円周率の近似値を求める
以下の図は,半角の公式を利用して,πの近似値を求めた例である.（πと近似値の値との絶対誤 差を表示している）
1e-09 1e-08 1e-07 1e-06 1e-05 0.0001 0.001 0.01 0.1 1
1 10 100 1000 10000 100000 1e+06 1e+07 1e+08 1e+09
Calculate Pi (1)
inner outer(1) outer(2)
/*
* 円周率の近似値を正多角形近似で求める
*/
#include <stdio.h>
#include <math.h>
#define MAX (24)
int main(int argc, char **argv) {
double inner, outer ; int i, n ;
n = 6 ; inner = 3.0 ;
outer = 2.0*sqrt(3.0) ;
printf("%10d %.15e %.15e\n", n, M_PI-inner, outer-M_PI) ; for(i=0;i<MAX;i++) {
/* ここで, 次のステップを計算する */
printf("%10d %.15e %.15e\n", n, fabs(M_PI-inner), fabs(outer-M_PI)) ; }
return 0 ; }
● 実習内容
（以下で「★」の数は推奨の程度を示します. 多いほど推奨の度合いが大きい. また,「†」は 難易度またはマニアックな程度を表します. 多いほど難しいかマニアックかです.）
1. （★★★）（前回の「4, 5」と同一）前回のプログラム１では, (1 + 1/n)n の計算を行うため に,実際にn−1回の乗算を行っていた. それを「繰り返し２乗法」（「高速乗算法」）および, 数学関数powを用いて書き直しなさい. また,絶対誤差のグラフを作成して,その誤差がn 回の乗算の場合とあまり差がないことを確かめなさい.
2. （★★★）半角の公式
ℓ2n= 2n s
1−p
1−(ℓn/n)2
2 ,
L2n= 2n2 Ln
(p
1 + (Ln/n)2−1)
を用いてπ の近似値を計算しなさい. さらに,πとそれぞれの近似値の絶対誤差のグラフを 書きなさい.
3. （★★★）半角の公式とそれから派生する公式
ℓ2n= 2n s
1−p
1−(ℓn/n)2
2 ,
2 L2n
= 1 ℓn + 1
Ln
を用いてπ の近似値を計算しなさい. さらに,πとそれぞれの近似値の絶対誤差のグラフを 書きなさい.
