余談ですが無理数の計算って何かワクワクしますよね。未知の数字を求めてる感じが凄く。
1.手計算で求める(その1)
一番簡単だと思う方法です。ただしあまり精度は良くないです。
まず2の10乗を計算してみましょう。
2の10乗までは大学入試で覚えておいたほうがいいと思います。でも計算すればすぐ出てきます。ここでこんな近似をしてします。
「≒」は大体等しいという意味です。1024は大体1000だろう,と大まかに決めてしまうのです。ここで両辺を底10の対数を取ります。
するとlog[10]2≒0.3という結果が出てきます。電卓で計算すると0.3010…なので小数点以下2桁が正しいことがわかります。
では他の場合はどうでしょう。3について考えてみます。さっきと同じように累乗していきます。
3^4=81はおおよそ80,つまり80=2^3✕10と表せます(3^2=9≒10じゃないかと思う人もいるかもしれませんが累乗した数は大きい方が精度が上がります)。両辺の対数をとると,
となって先ほど求めたlog[10]2 ≒ 0.3を代入すると
となります。log[10]3 = 0.47712なので小数点以下2桁が正しいことがわかります。
同様にして4,5,6,7,8,9,10,11,12…と求めていくことができますが,大きい数字になるに連れて誤差が大きくなっていくので気をつけてください。
基本的には素数さえ求めれば合成数は素因数の対数の和なので簡単に求められます。以下に5,7,11,13の場合を示します。
5の場合は10=5✕2ということを利用します(log[10]5=0.6989…)。
7の場合は7^2 = 49≒50を利用します(log[10]7=0.8450…)。
11の場合は11✕9 = 99 ≒ 100 を利用します。累乗してないので誤差は他に比べて大きいと思います(log[10]11=1.0413…)。
13の場合は13^3 = 2197 ≒ 2200 を利用します。この辺になるとどこで近似するかは好みですかね(log[10]13=1.1139…)。
極めると面白いものが見えてきそうですね。
2.手計算で求める(その2)
この方法はあまりおすすめはしません……が,こういう方法もあるんだと思ってみてってください。
対数は次の式によって表すことができます。ここで底はe(自然対数)です。
この式はf(x) = a^xを定義によって微分したものと,両辺の対数をとって微分した結果を照らし合わせると見えてきますが,循環論法になりそうなのでこれが成り立つものとして見て下さい。
ここで対数の重要な公式を使います。
(2)にa=10,c=e(自然対数),b=aを代入すると(1)より次の式が求まります。
a = 2について考えてみようと思うのですが,残念ながら2^h,10^hを求めるのがなかなか面倒なのです。h=0.5とすると,ルート2,ルート10の値を求めればよいのですが,
となり,正直うまくいきません。平方根ならばニュートン法や開平法で求められることには求められるのですが,正直あまり合理的ではなさそうですね。ちなみにh=1/16の場合で0.2860...です。
3.積分を使った方法
こっちは手計算で求めるというにはちょっと違うけれど,求められる方法。次の積分を用いる。
手計算で求めるならばこれを区分求積法で求めればオッケーです。でもごめんなさい。数学好きの高校生の子が見ていたら申し訳ないけれどここからプログラミングの話にシフトチェンジしていくのです。
上の式は積分ですが,積分を数値的に求める方法としてモンテカルロ法や台形公式なんかがまぁ楽でいいですね。それでちょっと実装してみました。ここでは自然対数で常用対数ではないのでそこんところよろしくお願いします。
/**
* 2015/10/13
* 対数を求めるプログラム
*/
using System;
class Program
{
// 試行回数
const int N = 10000;
static void Main()
{
for(int i = 2 ; i < 100 ; i++ )
{
Console.Write("log({0})\t",i);
Console.Write("{0:f5}\t",Math.Log(i));
Console.Write("{0:f5}\t",MyLog(i));
Console.Write("{0:f5}\t",MyLog2(i));
Console.WriteLine();
}
}
// モンテカルロ法で対数を求めるプログラム
static double MyLog(double d)
{
// 乱数生成用クラス
var rand = new Random();
var count = 0;
// 1未満の場合は定義しない
if( d < 1 )
return -1;
// 打点をする
for(int i = 0 ; i < N ; i++ )
{
// x = [1,d)の乱数生成
var x = rand.NextDouble() * ( d - 1 ) + 1;
// y = [0,1)の乱数生成
var y = rand.NextDouble();
// グラフの下側ならtrue
if( y < 1.0/x )
count++;
}
return (d-1.0) * count / N;
}
// 台形公式で対数を求めるプログラム
static double MyLog2(double d)
{
double integral = 0.0;
double dx = ( d - 1.0 )/N;
// 1未満の場合は定義しない
if( d < 1 )
return -1;
// 積分していく
for(int i = 1 ; i < N - 1 ; i++ )
integral += 1.0 / ( i * dx + 1);
// y0/2,yn/2を加える
integral += 1.0 / 2.0;
integral += 1.0 / d;
return integral * dx;
}
}
やっぱりモンテカルロ法には限界がありますね。
0 件のコメント:
コメントを投稿