直角三角形の斜辺を計算するアルゴリズム

もちろん、三平方の定理により平方根を求めれば良いだけなのですが、次のようなアルゴリズムにより四則演算のみで計算できます。Moler-Morrisonアルゴリズムというらしいですが、結果をみるとかなり収束が速いですね。

#include <stdio.h>
#include <stdlib.h>

/**
 * 直角三角形の斜辺の長さを計算する
 * @param a 直角を挟む辺1の長さ
 * @param b 直角を挟む辺2の長さ
 * @param iteration 繰り返し回数(増やすほど精度が上がる) 
 * @return 斜辺の長さ
 */
double mm_hypot(double a, double b, int iteration)
{
	double temp;
	int i;
	
	/* 絶対値 */
	a = (a>=0) ? a:-a;
	b = (b>=0) ? b:-b;
	
	/* 長辺をa, 短辺をb とする */
	if(a < b){
		temp = a;
		a = b;
		b = temp;
	}
	
	/* 短辺=0なら長辺を返す */
	if (b == 0) return a;

	/* 斜辺を求める繰り返し計算 (Moler-Morrisonアルゴリズム)*/
	for (i=0; i<iteration; i++){
		temp = (b/a) * (b/a);
		temp /= 4 + temp;
		a += 2*a*temp;
		b *= temp;
	}
	return a;
}

/**
 * テスト
 */
int main(void)
{
	double a,b,c;
	int i;

	a = 100;
	b = 100;
	printf("a=%f, b=%f\n",a,b);

	for(i=1;i<=5;i++){
		c = mm_hypot(a,b,i);
		printf("iteration = %d, c = %15.11f\n",i,c);
	}

	system("pause");
}



【結果】


a=100.000000, b=100.000000
iteration = 1, c = 140.00000000000
iteration = 2, c = 141.42131979695
iteration = 3, c = 141.42135623731
iteration = 4, c = 141.42135623731
iteration = 5, c = 141.42135623731