1 漸化式で定義される数列の計算

計算機数学演習 A 2015.5.19. Risa/Asir のプログラミング
前回 : Asir におけるベクトル (C 言語では配列に相当するもの)
漸化式で定義される数列の計算
1
例 1. 漸化式 an+1 = 3an + 2 (n ≥ 1), a1 = 1 で定義される数列 an (n = 1, 2, 3, . . .) を考える. この数列の第 1 項から
第 100 項までをベクトルに入れて, その値を表示するプログラムは次のようになる.
1. 漸化式で定義される数列の表示 (rec.txt)
A = newvect(101);
A[1] = 1;
for (I = 2; I <= 100; I++) {
A[I] = 3*A[I - 1] + 2;
}
print(A);
end;
漸化式 an+1 = 3an + 2 (n ≥ 1), a1 = 1 で定義される数列 an (n = 1, 2, 3, . . .) の第 1 項から第 N 項までをベクトル
に入れて, その値を表示する関数 rec は次のようになる. ここで N は関数 rec の引数である.
2. 漸化式で定義される数列を計算する関数 (rec2.txt)
def rec(N)
{
A = newvect(N + 1);
A[1] = 1;
for (I = 2; I <= N; I++) {
A[I] = 3*A[I - 1] + 2;
}
return A;
}
end;
問題 2. フィボナッチ数列とは, 漸化式 Fn+2 = Fn+1 + Fn (n ≥ 1), F1 = 1, F2 = 1 で定義される数列であった. フィボ
ナッチ数列 Fn の第 1 項から第 100 項までをベクトルに入れて, その値を表示するプログラムを書け.
3. フィボナッチ数列を表示 (fib.txt)
end;
問題 3. フィボナッチ数列 Fn の第 1 項から第 N 項までをベクトルに入れて, その値を表示する関数 fib を書け. ここ
で N は関数 fib の引数である.
4. フィボナッチ数列を計算する関数 (fib2.txt)
def fib(N)
{
}
end;
問題 4. チェビシェフ多項式 Tn (x) とは, 漸化式 Tn+1 (x) = 2xTn (x) − Tn−1 (x) (n ≥ 1), T0 (x) = 1, T1 (x) = x で定義
される多項式である. この多項式は cos の n 倍角の公式と関係があり, cos(nθ) = Tn (cos θ) が成り立つ.
1. チェビシェフ多項式 T2 (x), T3 (x), T4 (x), T5 (x) を計算せよ.
1
T2 (x)
T3 (x)
T4 (x)
T5 (x)
2. チェビシェフ多項式 T0 (x), T1 (x), T2 (x), . . . , T20 (x) をベクトルに入れて, その値を表示するプログラムを書け.
4. チェビシェフ多項式を表示 (cheb.txt)
end;
マクローリン展開の計算
2
関数 f (x) が x = 0 で n 階微分可能で連続であるとする. f (x) をマクローリン展開 (x = 0 のまわりのテイラー展開)
を行うと,
f (x) = f (0) + f ′ (0)x +
f (k) (0) k
f ′′ (0) 2
x + ··· +
x + ···
2!
k!
が得られる. 指数関数 ex のマクローリン展開は,
ex =
であった. Asir では無限和は扱えないので, 上のべき級数を x を N 次で打ち切った式を計算するプログラムを作る.
5. 指数関数 ex のマクローリン展開 (exp series.txt)
def exp_series(N)
{
S = 0;
for(K = 0; K <= N; K++) {
S = S + 1/fac(K) * x^K;
}
return S;
}
end;
/* fac(K) は階乗 K! を計算する関数 */
6. 指数関数 ex のマクローリン展開 (Asir 上)
(exp_series.txt を読み込む)
[224] F=exp_series(10);
<-- 変数 F に e^x を 10 次まで展開した多項式を代入
1/3628800*x^10+1/362880*x^9+1/40320*x^8+1/5040*x^7+1/720*x^6+1/120*x^5+1/24*x^4+1/6*x^3+1/2*x^2+x+1
[226] T=subst(F, x, 1);
<-- F の x に 1 を代入. (subst は式に値を代入するという関数) その結果を変数 T に代入.
9864101/3628800
[227] T*1.0;
<-- 分数を小数の値に直したい時は, 小数 1.0 をかけてやる. これは自然対数の底 e の値の近似値
2.71828
(もっと小数の精度を上げて計算したい場合)
[247] ctrl("bigfloat", 1);
<-- おまじない (大きい桁の小数を取り扱えるようにする)
[248] setprec(100);
<-- 小数 100 桁で計算させる
[249] F=exp_series(100)$
[250] T=subst(F, x, 1)$
[251] T*1.0;
2.7182818284590452353602874713526624977572470936999595749669676277240766303535475945713821785251664274274662
問題 5. tan x の逆関数 arctan x のマクローリン展開は
arctan x =
であるが, が, 上のべき級数を 2N + 1 次で打ち切った式を計算する関数 arctan_series(N) をつくれ. arctan(1) =
π
4
より円周率 π の近似値を計算せよ. また
1
π
1
− arctan
=
5
239
4
の公式 (マチンの公式, 1706 年, マチンはこれを使って円周率を 100 桁求めた) を用いて円周率 π の近似値を計算せよ.
4 arctan
2

Download Report