このアルゴリズムはJack Walderによって発明されました。JackWalderは、Convairで前述の爆撃機のナビゲーションコンピューターで働いていました。
桁ごとの方法の主な利点は、2つの操作による加算と除算のみを使用することです(右にシフトすることで簡単に実装できます)。
さらに、アルゴリズムは、ほとんどの計算機で使用されているバイナリ10進コードで直接動作するように作成できますが、以下の例では、このジャングルを使用しません。
このアルゴリズムは反復的であり、反復ごとに1つのアークタンジェントテーブルを使用します。 テーブルは事前に計算する必要があります。
#include <stdio.h> #include <math.h> int main(int argc, char **argv) { int bits = 32; int cordic_one = 1 << (bits - 2); printf("// , \n"); printf("static const int cordic_one = 0x%08x;\n", cordic_one); printf("static const int cordic_table[] = {\n"); double k = 1; for (int i = 0; i < bits; i++) { printf("0x%08x, // 0x%08x * atan(1/%.0f) \n", (int)(atan(pow(2, -i)) * cordic_one), cordic_one, pow(2, i)); k /= sqrt(1 + pow(2, -2 * i)); } printf("};\n"); printf("static const int cordic_k = 0x%08x; // %.16f * 0x%08x\n", (int)(k * cordic_one), k, cordic_one); }
同時に、スケーリング係数
cordic_k
が考慮されます。
その後、悪名高いsin 1°は次のように計算できます。
#include <stdio.h> #include <math.h> // , static const int cordic_one = 0x40000000; static const int cordic_table[] = { 0x3243f6a8, // 0x40000000 * atan(1/1) 0x1dac6705, // 0x40000000 * atan(1/2) 0x0fadbafc, // 0x40000000 * atan(1/4) 0x07f56ea6, // 0x40000000 * atan(1/8) 0x03feab76, // 0x40000000 * atan(1/16) 0x01ffd55b, // 0x40000000 * atan(1/32) 0x00fffaaa, // 0x40000000 * atan(1/64) 0x007fff55, // 0x40000000 * atan(1/128) 0x003fffea, // 0x40000000 * atan(1/256) 0x001ffffd, // 0x40000000 * atan(1/512) 0x000fffff, // 0x40000000 * atan(1/1024) 0x0007ffff, // 0x40000000 * atan(1/2048) 0x0003ffff, // 0x40000000 * atan(1/4096) 0x0001ffff, // 0x40000000 * atan(1/8192) 0x0000ffff, // 0x40000000 * atan(1/16384) 0x00007fff, // 0x40000000 * atan(1/32768) 0x00003fff, // 0x40000000 * atan(1/65536) 0x00001fff, // 0x40000000 * atan(1/131072) 0x00000fff, // 0x40000000 * atan(1/262144) 0x000007ff, // 0x40000000 * atan(1/524288) 0x000003ff, // 0x40000000 * atan(1/1048576) 0x000001ff, // 0x40000000 * atan(1/2097152) 0x000000ff, // 0x40000000 * atan(1/4194304) 0x0000007f, // 0x40000000 * atan(1/8388608) 0x0000003f, // 0x40000000 * atan(1/16777216) 0x0000001f, // 0x40000000 * atan(1/33554432) 0x0000000f, // 0x40000000 * atan(1/67108864) 0x00000008, // 0x40000000 * atan(1/134217728) 0x00000004, // 0x40000000 * atan(1/268435456) 0x00000002, // 0x40000000 * atan(1/536870912) 0x00000001, // 0x40000000 * atan(1/1073741824) 0x00000000, // 0x40000000 * atan(1/2147483648) }; static const int cordic_k = 0x26dd3b6a; // 0.6072529350088813 * 0x40000000 void cordic(int theta, int& s, int& c) { c = cordic_k; s = 0; for (int k = 0; k < 32; ++k) { int d = (theta >= 0) ? 0 : -1; int tx = c - (((s >> k) ^ d) - d); int ty = s + (((c >> k) ^ d) - d); c = tx; s = ty; theta -= ((cordic_table[k] ^ d) - d); } } int main(void) { double alpha = M_PI / 180; int sine, cosine; cordic(alpha * cordic_one, sine, cosine); printf("CORDIC: %.8f\nExpected: %.8f\n", (double)sine / cordic_one, sin(alpha)); }
結果:
CORDIC: 0.01745240 Expected: 0.01745241
32回の反復があるため、小さなエラーが残ります。 計算機は通常40回の反復を使用します。