べき乗の計算について

質問

ｘ＾０．５　　をやりたいのですが　＾が実施できないコンパイラを使っています。
ｘ＾３　なら　ｘ＊ｘ＊ｘとすればいいのですが　＾０．５の近似式を教えて
欲しいのですが・・

回答No.4

stomachman · ベストアンサー

√x = x^0.5　は言語によってはsqrt(x), sqr(x), power(x,0.5) などと書かれます。

そういうのがない場合には、数値計算をする関数（functionあるいはprocedure, subroutine, method）を作ってやれば良いですね。 要するに、大抵の言語でコンパイラに作りつけの関数sqrtがやっていることを自分で書けばよい。

～～～～～～～～～～～～～～～～～～～～～

√xを高い精度で計算するには、漸化式を使うのが速いと思います。

もちろんx<0のときはエラー、x=0なら√x=0です。
x>0のとき、y=√x に近づく近似値の列y[0], y[1], ...., y[k], .... を次々と作り、誤差が小さくなったら終わります。

●漸化式
適当な出発値y[0]から始めて
y[n+1]＝y[n]/2+x/(2y[n])
という漸化式でy[1], y[2], ... を作ります。（Newton法）
この式は、繰り返しの度に有効桁数が倍になる、「２次収束」という性質を持っているので、出発値y[0]が正解に近ければ、ほんの数回の繰り返しで十分な精度（10進法で10桁以上）に達します。

●繰り返しを何処で打ち切るかの判断
　筋から言うと、正解√xに対する相対誤差 δ
δ=｜(y[n]^2)/x - 1｜/2
が十分小さくなったら打ち切る、ということで良いでしょう。
　このδがなぜ相対誤差なのか。近似値y[n]が絶対誤差εを含んでいて
y[n]=sqrt(x)+ε
であるとすると、
δ=｜(y[n]^2)/x - 1｜/2
= ｜(x+2y[n]ε+(ε^2))/x - 1｜/2
= ｜(2ε√x+ε^2)/x ｜/2
≒｜(2ε√x)/x ｜/2　　　（|ε|は小さいのでε^2 は無視できる。)
= ｜(ε/√x ｜
だからです。
　幾らなら十分小さくなったと言えるか。これは要求される精度に依存しますが、一般に計算機内部での浮動小数点数値の有効桁数を超えたら、それ以上計算したってしょうがないですね。
　これを逆手に取って、
　　y[n]とy[n+1]が同じになったら打ち切る
というやり方もあります。

●出発値y[0]
出発値y[0]があまりにもでたらめだと、上記の漸化式が旨く働かないことがあります。
　y[0]は、数値xの計算機内部での表現を利用すると簡単に決められます。
浮動小数点の数値xは符号部S、仮数部D、指数部Eで表され、たいてい１６進法で
x = S×D×(16^E)
となっている。Eは整数です。
ここで、x≠0の場合には
1＞D≧1/16
です。
　この問題ではx≠0だし、符号部SはS=1に決まっています。従って
√x = √D ×(4^E)
ですね。
√D=√(1+(D-1))≒1+(D-1)/2　≒1
と荒っぽい近似ができます。また
Eが偶数の場合には
√x　= √D ×16^(E/2)
Eが奇数の場合には
√x = (√D )/4 ×16^((E+1)/2)
です。

だから
　　Eが偶数のとき　y[0] = 符号部1 仮数部(1+(D-1)/2) 指数部(E/2)
　　Eが奇数のとき　y[0] = 符号部1 仮数部(1+(D-1)/2)/4 指数部((E+1)/2)
とでもすれば良いでしょう。もっと手抜きして
　　Eが偶数のとき　y[0] = 符号部1 仮数部1 指数部(E/2)
　　Eが奇数のとき　y[0] = 符号部1 仮数部1 指数部((E+1)/2)
でも大丈夫です。（漸化式の繰り返しが少し増えるだけです。）

内部表現の調べかたが分からない場合には、計算でDとEを出すこともできます。
　16^(E-1) < x < 16^(E)
となるEを見つければ良いのです。

お礼

ごていねいに　ありがとうございました。
結局、ｓｑｒｔ関数を手作りでつくることにしました。

回答No.3

guiter

sight さんが仰るように何か関数があるような気はしますが、
どうしても近似が必要であれば次のテイラー展開の結果を使ってください。

　(1+x)^α = 1 + (α,1)*x + (α,2)*x^2 + … + (α,n)*x^n + o(x^n)

ここで、(α,n) は２項係数でαは自然数でなくてもよいですが n は自然数です。
（通常カッコの中のαと n は縦に並べて書きますが今は横に並べています。）
つまり、

　(α,n) = α(α-1)…(α-n+1)/n!

です。ただし、(α,0)=1 です。
今は、α=1/2 ですから
　(1/2,1) = (1/2)/1! = 1/2
　(1/2,2) = (1/2)(1/2-1)/2! = -1/8
　(1/2,3) = (1/2)(1/2-1)(1/2-2)/3! = 1/16
などになります。
したがって、

　(1+x)^0.5 = 1 + x/2 - x^2 /8 + x^3 /16 …

と展開できます。ご質問の x^0.5 であれば、x+1 ⇒ x とずらして

　x^0.5 = 1 + (x-1)/2 - (x-1)^2 /8 + (x-1)^3 /16 …

ですね。
あとは近似の精度に気をつけて下さい。

回答No.2

sight

コンパイラ・・・プログラムですか？
^が使えなくても、Cのpow関数みたいに、べき剰を計算してくれる関数はないんでしょうか？
pow(x, 0.5)みたいにかけるような・・・・

回答No.1

k_eba

．５乗だけでよいなら
ルートを使ってください。

また

参考で色々書いていますよ

参考URL：

http://oshiete1.goo.ne.jp/kotaeru.php3?q=144106