CORDIC: Difference between revisions

31,893 bytes added ,  1 month ago
m
→‎{{header|Phix}}: use pygments
m (→‎Notes: tighten up estimate on observed precision problems)
m (→‎{{header|Phix}}: use pygments)
 
(16 intermediate revisions by 10 users not shown)
Line 1:
{{draftDraft task|Trigonometrics}}
 
;Introduction
Line 41:
 
REAL pi by 2 = pi / 2;
REAL pi by 4 = pi / 4;
 
# pre-computed table of arctan(10^-n) values: #
[]REAL theta = []REAL( 7.85398163397448e-01, 9.96686524911620e-02
, 9.99966668666524e-03, 9.99999666666867e-04
, 9.99999996666667e-05, 9.99999999966667e-06
, 9.99999999999667e-07, 9.99999999999997e-08
, 1.00000000000000e-08, 1.00000000000000e-09
, 1.00000000000000e-10, 1.00000000000000e-11
, 1.00000000000000e-12, 1.00000000000000e-13
, 1.00000000000000e-14, )[ AT 0 ];1.00000000000000e-15
, 1.00000000000000e-16
REAL epsilon = 1e-8;
);
REAL epsilon = 1e-16;
 
# mode to hold the values returned from the CORDIC procedure #
Line 83 ⟶ 84:
IF flip sign AND sign both THEN sign := -sign FI;
REAL x := 1, y := 0;
INT k := 01;
REAL ten to minus k := 1; # NB: ten to minus k #= 10.0^0-(k-1) #
WHILE epsilon < alpha DO
WHILE alpha < theta[ k ] DO
Line 128 ⟶ 129:
REAL s diff = ABS ( cordic sine - sine );
REAL t diff = ABS ( cordic tan - tangent );
print( ( fixed( angle, -48, 14 ), ": "
, fixed( cordic cosine, -9, 6 ), " "
, fixed( cosine, -9, 6 ), " "
Line 143 ⟶ 144:
 
[]REAL tests = ( -9, 0, 1.5, 6 );
print( ( " angle cordic cos cos difference" ) );
print( ( " cordic sin sin difference" ) );
print( ( " cordic tan tan difference" ) );
print( ( newline ) );
FOR i FROM LWB tests TO UPB tests DO
show cordic( tests[ i ] )
OD;
print( ( " angle cordic cos cos difference" ) );
print( ( " cordic sin sin difference" ) );
print( ( " cordic tan tan difference" ) );
print( ( newline ) );
FOR i FROM 7810 BY 9 TO 7898 DO
show cordic( i / 10 000 )
OD
 
Line 155 ⟶ 163:
{{out}}
<pre>
angle cordic cos cos difference cordic sin sin difference cordic tan tan difference
-9.00000: -0.911130 -0.911130 31.88699162e1102230e-916 | -0.412118 -0.412118 80.59353766e-900000000e+0 | 0.452316 0.452316 10.13613664e-800000000e+0
0.00000: 1.000000 1.000000 0.00000000e+0 | 0.000000 0.000000 0.00000000e+0 | 0.000000 0.000000 0.00000000e+0
1.55000: 0.070737 0.070737 47.46952739e4940054e-916 | 0.997495 0.997495 30.1695591e-1000000000e+0 | 14.101419101420 14.101420 81.95478376e5099033e-713
6.00000: 0.960170 0.960170 21.71217671e1102230e-916 | -0.279415 -0.279415 95.31999689e5511151e-916 | -0.291006 -0.291006 16.05286085e1062266e-816
angle cordic cos cos difference cordic sin sin difference cordic tan tan difference
0.7810: 0.710210 0.710210 3.3306691e-16 | 0.703990 0.703990 2.2204460e-16 | 0.991242 0.991242 6.6613381e-16
0.7819: 0.709576 0.709576 2.2204460e-16 | 0.704629 0.704629 1.1102230e-16 | 0.993028 0.993028 4.4408921e-16
0.7828: 0.708942 0.708942 3.3306691e-16 | 0.705267 0.705267 4.4408921e-16 | 0.994817 0.994817 9.9920072e-16
0.7837: 0.708307 0.708307 6.6613381e-16 | 0.705905 0.705905 5.5511151e-16 | 0.996609 0.996609 1.6653345e-15
0.7846: 0.707671 0.707671 5.5511151e-16 | 0.706542 0.706542 3.3306691e-16 | 0.998405 0.998405 1.3322676e-15
0.7855: 0.707035 0.707035 1.1102230e-16 | 0.707179 0.707179 1.1102230e-16 | 1.000204 1.000204 2.2204460e-16
0.7864: 0.706398 0.706398 1.1102230e-16 | 0.707815 0.707815 2.2204460e-16 | 1.002006 1.002006 4.4408921e-16
0.7873: 0.705761 0.705761 1.1102230e-16 | 0.708450 0.708450 1.1102230e-16 | 1.003811 1.003811 2.2204460e-16
0.7882: 0.705123 0.705123 3.3306691e-16 | 0.709085 0.709085 5.5511151e-16 | 1.005619 1.005619 1.3322676e-15
0.7891: 0.704484 0.704484 5.5511151e-16 | 0.709720 0.709720 4.4408921e-16 | 1.007431 1.007431 1.3322676e-15
</pre>
 
Line 412 ⟶ 431:
 
Also: [[j:Vocabulary/SpecialCombinations#Preexecuting_verbs_with_.28.28_.29.29|double parenthesis around a noun phrase]] tells the interpreter that that expression is a constant which should be evaluated once, ahead of time.
 
=={{header|Java}}==
{{trans|C}}
<syntaxhighlight lang="Java">
public class CordicCalculator {
 
// Constants
private static final double[] angles = {
0.78539816339745, 0.46364760900081, 0.24497866312686, 0.12435499454676,
0.06241880999596, 0.03123983343027, 0.01562372862048, 0.00781234106010,
0.00390623013197, 0.00195312251648, 0.00097656218956, 0.00048828121119,
0.00024414062015, 0.00012207031189, 0.00006103515617, 0.00003051757812,
0.00001525878906, 0.00000762939453, 0.00000381469727, 0.00000190734863,
0.00000095367432, 0.00000047683716, 0.00000023841858, 0.00000011920929,
0.00000005960464, 0.00000002980232, 0.00000001490116, 0.00000000745058
};
 
private static final double[] kvalues = {
0.70710678118655, 0.63245553203368, 0.61357199107790, 0.60883391251775,
0.60764825625617, 0.60735177014130, 0.60727764409353, 0.60725911229889,
0.60725447933256, 0.60725332108988, 0.60725303152913, 0.60725295913894,
0.60725294104140, 0.60725293651701, 0.60725293538591, 0.60725293510314,
0.60725293503245, 0.60725293501477, 0.60725293501035, 0.60725293500925,
0.60725293500897, 0.60725293500890, 0.60725293500889, 0.60725293500888
};
 
// Convert degrees to radians
private static double radians(double degrees) {
return degrees * Math.PI / 180.0;
}
 
// Cordic algorithm implementation
private static void cordic(double alpha, int n, double[] result) {
int i, ix, sigma;
double kn, x, y, atn, t, theta = 0.0, pow2 = 1.0;
int newsgn = (int)Math.floor(alpha / (2.0 * Math.PI)) % 2 == 1 ? 1 : -1;
if (alpha < -Math.PI/2.0 || alpha > Math.PI/2.0) {
if (alpha < 0) {
cordic(alpha + Math.PI, n, result);
} else {
cordic(alpha - Math.PI, n, result);
}
result[0] = result[0] * newsgn;
result[1] = result[1] * newsgn;
return;
}
ix = n - 1;
if (ix > 23) ix = 23;
kn = kvalues[ix];
x = 1;
y = 0;
for (i = 0; i < n; ++i) {
atn = angles[i];
sigma = (theta < alpha) ? 1 : -1;
theta += sigma * atn;
t = x;
x -= sigma * y * pow2;
y += sigma * t * pow2;
pow2 /= 2.0;
}
result[0] = x * kn;
result[1] = y * kn;
}
 
// Main method
public static void main(String[] args) {
double c_cos, c_sin;
double[] result = new double[2];
double[] angles = {-9.0, 0.0, 1.5, 6.0};
String format;
System.out.println(" x sin(x) diff. sine cos(x) diff. cosine");
format = "%+03d.0° %+.8f (%+.8f) %+.8f (%+.8f)\n";
for (int th = -90; th <= +90; th += 15) {
double thr = radians(th);
cordic(thr, 24, result);
c_cos = result[0];
c_sin = result[1];
System.out.printf(format, th, c_sin, c_sin - Math.sin(thr), c_cos, c_cos - Math.cos(thr));
}
 
System.out.println("\nx(rads) sin(x) diff. sine cos(x) diff. cosine");
format = "%+4.1f %+.8f (%+.8f) %+.8f (%+.8f)\n";
for (int i = 0; i < angles.length; ++i) {
double thr = angles[i];
cordic(thr, 24, result);
c_cos = result[0];
c_sin = result[1];
System.out.printf(format, thr, c_sin, c_sin - Math.sin(thr), c_cos, c_cos - Math.cos(thr));
}
}
}
</syntaxhighlight>
{{out}}
<pre>
x sin(x) diff. sine cos(x) diff. cosine
-90.0° -1.00000000 (+0.00000000) -0.00000007 (-0.00000007)
-75.0° -0.96592585 (-0.00000003) +0.25881895 (-0.00000009)
-60.0° -0.86602545 (-0.00000005) +0.49999992 (-0.00000008)
-45.0° -0.70710684 (-0.00000006) +0.70710672 (-0.00000006)
-30.0° -0.49999992 (+0.00000008) +0.86602545 (+0.00000005)
-15.0° -0.25881895 (+0.00000009) +0.96592585 (+0.00000003)
+00.0° +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+15.0° +0.25881895 (-0.00000009) +0.96592585 (+0.00000003)
+30.0° +0.49999992 (-0.00000008) +0.86602545 (+0.00000005)
+45.0° +0.70710684 (+0.00000006) +0.70710672 (-0.00000006)
+60.0° +0.86602545 (+0.00000005) +0.49999992 (-0.00000008)
+75.0° +0.96592585 (+0.00000003) +0.25881895 (-0.00000009)
+90.0° +1.00000000 (-0.00000000) -0.00000007 (-0.00000007)
 
x(rads) sin(x) diff. sine cos(x) diff. cosine
-9.0 -0.41211842 (+0.00000006) -0.91113029 (-0.00000003)
+0.0 +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+1.5 +0.99749499 (+0.00000000) +0.07073719 (-0.00000002)
+6.0 -0.27941552 (-0.00000002) +0.96017028 (-0.00000001)
 
</pre>
 
 
=={{header|jq}}==
'''Adapted from [[#Wren|Wren]]'''
{{works with|jq}}
 
'''Works with gojq, the Go implementation of jq'''
 
'''Preliminaries'''
<syntaxhighlight lang="jq">
# like while/2 but emit the final term rather than the first one
def whilst(cond; update):
def _whilst:
if cond then update | (., _whilst) else empty end;
_whilst;
 
# Simplistic approach to rounding:
def round($ndec): pow(10;$ndec) as $p | . * $p | round / $p;
 
def lpad($len): tostring | ($len - length) as $l | (" " * $l) + .;
</syntaxhighlight>
'''The task'''
<syntaxhighlight lang="jq">
# The following are pre-computed to avoid using atan and sqrt functions.
def angles: [
0.78539816339745, 0.46364760900081, 0.24497866312686, 0.12435499454676,
0.06241880999596, 0.03123983343027, 0.01562372862048, 0.00781234106010,
0.00390623013197, 0.00195312251648, 0.00097656218956, 0.00048828121119,
0.00024414062015, 0.00012207031189, 0.00006103515617, 0.00003051757812,
0.00001525878906, 0.00000762939453, 0.00000381469727, 0.00000190734863,
0.00000095367432, 0.00000047683716, 0.00000023841858, 0.00000011920929,
0.00000005960464, 0.00000002980232, 0.00000001490116, 0.00000000745058
];
 
def kvalues: [
0.70710678118655, 0.63245553203368, 0.61357199107790, 0.60883391251775,
0.60764825625617, 0.60735177014130, 0.60727764409353, 0.60725911229889,
0.60725447933256, 0.60725332108988, 0.60725303152913, 0.60725295913894,
0.60725294104140, 0.60725293651701, 0.60725293538591, 0.60725293510314,
0.60725293503245, 0.60725293501477, 0.60725293501035, 0.60725293500925,
0.60725293500897, 0.60725293500890, 0.60725293500889, 0.60725293500888
];
 
def PI: (1 | atan * 4);
 
def radians: . * PI / 180;
 
def Cordic($alpha; $n):
PI as $PI
| (if (($alpha / (2 * $PI))|floor % 2 == 1) then 1 else -1 end) as $newsgn
| if ($alpha < -$PI/2)
then Cordic($alpha + $PI; $n) | map(. * $newsgn)
elif ($alpha > $PI/2)
then Cordic($alpha - $PI; $n) | map(. * $newsgn)
else kvalues[ [$n-1, (kvalues|length-1)] | min] as $kn
| reduce angles[0:$n][] as $atan ({ theta: 0, x: 1, y: 0, pow2: 1};
(if .theta < $alpha then 1 else -1 end) as $sigma
| .theta += $sigma * $atan
| .x as $t
| .x = (.x - $sigma * .y * .pow2)
| .y = (.y + $sigma * $t * .pow2)
| .pow2 /= 2 )
| [.x * $kn, .y * $kn]
end;
 
def example:
def format: map(round(8) | lpad(11)) | join(" ");
" x sin(x) diff. sine cos(x) diff. cosine",
({ th: -90 }
| whilst (.th <= 90;
(.th | radians) as $thr
| Cordic($thr; 24) as [ $cos, $sin ]
| .out = [.th, $sin, $sin - ($thr|sin), $cos, $cos - ($thr|cos) ]
| .th += 15 )
| .out
| format),
 
"\n x(rads) sin(x) diff. sine cos(x) diff. cosine",
((-9, 0, 1.5, 6) as $thr
| Cordic($thr; 24) as [$cos, $sin]
| [$thr, $sin, $sin - ($thr|sin), $cos, $cos - ($thr|cos)]
| format) ;
 
example
</syntaxhighlight>
'''Invocation:''' jq -ncr -f cordic.jq
{{output}}
<pre>
x sin(x) diff. sine cos(x) diff. cosine
-90 -1 0 -7e-08 -7e-08
-75 -0.96592585 -3e-08 0.25881895 -9e-08
-60 -0.86602545 -5e-08 0.49999992 -8e-08
-45 -0.70710684 -6e-08 0.70710672 -6e-08
-30 -0.49999992 8e-08 0.86602545 5e-08
-15 -0.25881895 9e-08 0.96592585 3e-08
0 7e-08 7e-08 1 -0
15 0.25881895 -9e-08 0.96592585 3e-08
30 0.49999992 -8e-08 0.86602545 5e-08
45 0.70710684 6e-08 0.70710672 -6e-08
60 0.86602545 5e-08 0.49999992 -8e-08
75 0.96592585 3e-08 0.25881895 -9e-08
90 1 -0 -7e-08 -7e-08
 
x(rads) sin(x) diff. sine cos(x) diff. cosine
-9 -0.41211842 6e-08 -0.91113029 -3e-08
0 7e-08 7e-08 1 -0
1.5 0.99749499 0 0.07073719 -2e-08
6 -0.27941552 -2e-08 0.96017028 -1e-08
</pre>
 
=={{header|Julia}}==
Line 522 ⟶ 767:
+6.0 -0.27941552 (-0.00000002) +0.96017028 (-0.00000001)
</pre>
 
=={{header|Lua}}==
{{Trans|ALGOL 68}}
<syntaxhighlight lang="lua">
do -- implement sin, cos and tan using the CORDIC algorithm
 
local piBy2 = math.pi / 2
 
-- pre-computed table of arctan(10^-n) values
local theta = { 7.85398163397448e-01, 9.96686524911620e-02
, 9.99966668666524e-03, 9.99999666666867e-04
, 9.99999996666667e-05, 9.99999999966667e-06
, 9.99999999999667e-07, 9.99999999999997e-08
, 1.00000000000000e-08, 1.00000000000000e-09
, 1.00000000000000e-10, 1.00000000000000e-11
, 1.00000000000000e-12, 1.00000000000000e-13
, 1.00000000000000e-14, 1.00000000000000e-15
, 1.00000000000000e-16
}
local epsilon = 1e-16
 
--[[ CORDIC algorithm, finds "y" and "x" for alpha radians
signs indicates the sign of the results:
signs[ 1 ] = sign for -π : -π/2, signs[ 1 ] = sign for -π/2 : 0
signs[ 3 ] = sign for 0 : π/2, signs[ 4 ] = sign for π/2 : π
signBoth = TRUE => sign applied to both y and x, FALSE => sign y only
--]]
function cordic( alphaIn, signs, signBoth )
local alpha = alphaIn
-- ensure -π <= alpha <= π
local flipSign = false
while alpha < - math.pi do alpha = alpha + math.pi flipSign = not flipSign end
while alpha > math.pi do alpha = alpha - math.pi flipSign = not flipSign end
local sign
if alpha < - piBy2 then
alpha = alpha + math.pi
sign = signs[ 1 ]
elseif alpha < 0 then
alpha = - alpha
sign = signs[ 2 ]
elseif alpha < piBy2 then
sign = signs[ 3 ]
else -- alpha <= math.pi
alpha = math.pi - alpha
sign = signs[ 4 ]
end
if flipSign and signBoth then sign = -sign end
local x, y, k, tenToMinusK = 1, 0, 1, 1 -- NB: tenToMinusK is 10^-(k-1)
while epsilon < alpha do
while alpha < theta[ k ] do
k = k + 1
tenToMinusK = tenToMinusK / 10
end
alpha = alpha - theta[ k ]
x, y = x - tenToMinusK * y, y + tenToMinusK * x
end
return sign * y, signBoth and sign * x or x
end
 
-- cos(alpha) using the CORDIC algorithm. alpha in radians
function cCos( alpha )
local y, x = cordic( alpha, { -1, 1, 1, -1 }, true )
return x / math.sqrt( ( x * x ) + ( y * y ) )
end
 
-- sin(alpha) using the CORDIC algorithm. alpha in radians
function cSin( alpha )
local y, x = cordic( alpha, { -1,-1, 1, 1 }, true )
return y / math.sqrt( ( x * x ) + ( y * y ) )
end
 
-- tan(alpha) using the CORDIC algorithm. alpha in radians
function cTan( alpha )
local y, x = cordic( alpha, { 1, -1, 1,-1 }, false )
return x == 0 and NaN or y / x
end
 
function f6( v ) return string.format( " %9.6f", v ) end
function g( v ) return string.format( " %12g", v ) end
 
function showCordic( angle )
local cosine, cordicCos = math.cos( angle ), cCos( angle )
local sine, cordicSin = math.sin( angle ), cSin( angle )
local tangent, cordicTan = math.tan( angle ), cTan( angle )
local cDiff = cordicCos - cosine if cDiff < 0 then cDiff = - cDiff end
local sDiff = cordicSin - sine if sDiff < 0 then sDiff = - sDiff end
local tDiff = cordicTan - tangent if tDiff < 0 then tDiff = - tDiff end
io.write( string.format( "%4.1f: ", angle )
, f6( cordicCos ), f6( cosine ), g( cDiff ), " |"
, f6( cordicSin ), f6( sine ), g( sDiff ), " |"
, f6( cordicTan ), f6( tangent ), g( tDiff ), "\n"
)
end
 
local tests = { -9, 0, 1.5, 6 }
io.write( "angle cordic cos cos difference" )
io.write( " cordic sin sin difference" )
io.write( " cordic tan tan difference" )
io.write( "\n" )
for i = 1, # tests do
showCordic( tests[ i ] )
end
io.write( "angle cordic cos cos difference" )
io.write( " cordic sin sin difference" )
io.write( " cordic tan tan difference" )
io.write( "\n" )
for i = 7810, 7891, 9 do
showCordic( i / 10000 )
end
 
end
</syntaxhighlight>
{{out}}
<pre>
angle cordic cos cos difference cordic sin sin difference cordic tan tan difference
-9.0: -0.911130 -0.911130 1.11022e-16 | -0.412118 -0.412118 0 | 0.452316 0.452316 0
0.0: 1.000000 1.000000 0 | 0.000000 0.000000 0 | 0.000000 0.000000 0
1.5: 0.070737 0.070737 7.49401e-16 | 0.997495 0.997495 0 | 14.101420 14.101420 1.5099e-13
6.0: 0.960170 0.960170 1.11022e-16 | -0.279415 -0.279415 5.55112e-16 | -0.291006 -0.291006 6.10623e-16
angle cordic cos cos difference cordic sin sin difference cordic tan tan difference
0.8: 0.710210 0.710210 3.33067e-16 | 0.703990 0.703990 2.22045e-16 | 0.991242 0.991242 6.66134e-16
0.8: 0.709576 0.709576 2.22045e-16 | 0.704629 0.704629 1.11022e-16 | 0.993028 0.993028 4.44089e-16
0.8: 0.708942 0.708942 3.33067e-16 | 0.705267 0.705267 4.44089e-16 | 0.994817 0.994817 9.99201e-16
0.8: 0.708307 0.708307 6.66134e-16 | 0.705905 0.705905 5.55112e-16 | 0.996609 0.996609 1.66533e-15
0.8: 0.707671 0.707671 5.55112e-16 | 0.706542 0.706542 3.33067e-16 | 0.998405 0.998405 1.33227e-15
0.8: 0.707035 0.707035 1.11022e-16 | 0.707179 0.707179 1.11022e-16 | 1.000204 1.000204 2.22045e-16
0.8: 0.706398 0.706398 1.11022e-16 | 0.707815 0.707815 2.22045e-16 | 1.002006 1.002006 4.44089e-16
0.8: 0.705761 0.705761 1.11022e-16 | 0.708450 0.708450 1.11022e-16 | 1.003811 1.003811 2.22045e-16
0.8: 0.705123 0.705123 3.33067e-16 | 0.709085 0.709085 5.55112e-16 | 1.005619 1.005619 1.33227e-15
0.8: 0.704484 0.704484 5.55112e-16 | 0.709720 0.709720 4.44089e-16 | 1.007431 1.007431 1.33227e-15
</pre>
 
=={{header|Perl}}==
{{trans|Raku}}
<syntaxhighlight lang="perl" line>
use strict;
use warnings;
 
sub CORDIC {
my($a) = shift;
my ($k, $x, $y) = (0, 1, 0);
my @PoT = (1, 0.1, 0.01, 0.001, 0.0001, 0.00001);
my @Tbl = <7.853981633e-1 9.966865249e-2 9.999666686e-3 9.999996666e-4 9.999999966e-5 9.999999999e-6 0>;
 
while ($a > 1e-5) {
$k++ while $a < $Tbl[$k];
$a -= $Tbl[$k];
($x,$y) = ($x - $PoT[$k]*$y, $y + $PoT[$k]*$x);
}
$x / sqrt($x*$x + $y*$y)
}
 
print "Angle CORDIC Cosine Error\n";
for my $angle (<-9 0 1.5 6>) {
my $cordic = CORDIC abs $angle;
printf "%4.1f %12.8f %12.8f %12.8f\n", $angle, $cordic, cos($angle), cos($angle) - $cordic
}
</syntaxhighlight>
{{out}}
<pre>
Angle CORDIC Cosine Error
-9.0 -0.91112769 -0.91113026 -0.00000257
0.0 1.00000000 1.00000000 0.00000000
1.5 0.07073880 0.07073720 -0.00000160
6.0 0.96016761 0.96017029 0.00000268
</pre>
 
=={{header|Phix}}==
{{trans|Wren}}
<!--(phixonline)-->
<syntaxhighlight lang="phix">
with javascript_semantics
constant angles = {
0.78539816339745, 0.46364760900081, 0.24497866312686, 0.12435499454676,
0.06241880999596, 0.03123983343027, 0.01562372862048, 0.00781234106010,
0.00390623013197, 0.00195312251648, 0.00097656218956, 0.00048828121119,
0.00024414062015, 0.00012207031189, 0.00006103515617, 0.00003051757812,
0.00001525878906, 0.00000762939453, 0.00000381469727, 0.00000190734863,
0.00000095367432, 0.00000047683716, 0.00000023841858, 0.00000011920929,
0.00000005960464, 0.00000002980232, 0.00000001490116, 0.00000000745058}
constant kvalues = {
0.70710678118655, 0.63245553203368, 0.61357199107790, 0.60883391251775,
0.60764825625617, 0.60735177014130, 0.60727764409353, 0.60725911229889,
0.60725447933256, 0.60725332108988, 0.60725303152913, 0.60725295913894,
0.60725294104140, 0.60725293651701, 0.60725293538591, 0.60725293510314,
0.60725293503245, 0.60725293501477, 0.60725293501035, 0.60725293500925,
0.60725293500897, 0.60725293500890, 0.60725293500889, 0.60725293500888}
 
function cordic(atom alpha, integer n)
bool sgn = +1
while alpha < -PI/2 do alpha += PI sgn *= -1 end while
while alpha > PI/2 do alpha -= PI sgn *= -1 end while
atom kn = sgn * kvalues[min(n,length(kvalues))],
atan, theta = 0, x = 1, y = 0, pow2 = 1
for i=1 to n do
atan = iff(i<=length(angles)?angles[i]:atan/2)
atom sigma = iff(theta < alpha ? 1 : -1), t = x
-- atom sigma = iff(theta <= alpha ? 1 : -1), t = x -- (matches Julia)
theta += sigma * atan
x -= sigma * y * pow2
y += sigma * t * pow2
pow2 /= 2
end for
return {x * kn, y * kn}
end function
 
constant fmh = "%n x(%s) sin(x) diff.sine cos(x) diff.cosine\n",
fmt = "%+6.1f %+.8f (%+.8f) %+.8f (%+.8f)\n"
printf(1,fmh,{false,"deg"})
for theta = -90 to 90 by 15 do
atom r = theta*PI/180,
{c,s} = cordic(r,24)
printf(1,fmt,{theta,s,s-sin(r),c,c-cos(r)})
end for
printf(1,fmh,{true,"rad"})
for r in {-9, 0, 1.5, 6} do
atom {c,s} = cordic(r, 24)
printf(1,fmt,{r,s,s-sin(r),c,c-cos(r)})
end for
</syntaxhighlight>
{{out}}
As noted above in the commented out line, and mentioned in the Wren entry, testing for <= flips the sign of sin(0).
<pre>
x(deg) sin(x) diff.sine cos(x) diff.cosine
-90.0 -1.00000000 (+0.00000000) -0.00000007 (-0.00000007)
-75.0 -0.96592585 (-0.00000003) +0.25881895 (-0.00000009)
-60.0 -0.86602545 (-0.00000005) +0.49999992 (-0.00000008)
-45.0 -0.70710684 (-0.00000006) +0.70710672 (-0.00000006)
-30.0 -0.49999992 (+0.00000008) +0.86602545 (+0.00000005)
-15.0 -0.25881895 (+0.00000009) +0.96592585 (+0.00000003)
+0.0 +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+15.0 +0.25881895 (-0.00000009) +0.96592585 (+0.00000003)
+30.0 +0.49999992 (-0.00000008) +0.86602545 (+0.00000005)
+45.0 +0.70710684 (+0.00000006) +0.70710672 (-0.00000006)
+60.0 +0.86602545 (+0.00000005) +0.49999992 (-0.00000008)
+75.0 +0.96592585 (+0.00000003) +0.25881895 (-0.00000009)
+90.0 +1.00000000 (-0.00000000) -0.00000007 (-0.00000007)
 
x(rad) sin(x) diff.sine cos(x) diff.cosine
-9.0 -0.41211842 (+0.00000006) -0.91113029 (-0.00000003)
+0.0 +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+1.5 +0.99749499 (+0.00000000) +0.07073719 (-0.00000002)
+6.0 -0.27941552 (-0.00000002) +0.96017028 (-0.00000001)
</pre>
 
=={{header|Python}}==
{{trans|C}}
<syntaxhighlight lang="python">
# CORDIC.py by Xing216
from math import pi, sin, cos, floor
 
angles = [
0.78539816339745, 0.46364760900081, 0.24497866312686, 0.12435499454676,
0.06241880999596, 0.03123983343027, 0.01562372862048, 0.00781234106010,
0.00390623013197, 0.00195312251648, 0.00097656218956, 0.00048828121119,
0.00024414062015, 0.00012207031189, 0.00006103515617, 0.00003051757812,
0.00001525878906, 0.00000762939453, 0.00000381469727, 0.00000190734863,
0.00000095367432, 0.00000047683716, 0.00000023841858, 0.00000011920929,
0.00000005960464, 0.00000002980232, 0.00000001490116, 0.00000000745058
]
kvalues = [
0.70710678118655, 0.63245553203368, 0.61357199107790, 0.60883391251775,
0.60764825625617, 0.60735177014130, 0.60727764409353, 0.60725911229889,
0.60725447933256, 0.60725332108988, 0.60725303152913, 0.60725295913894,
0.60725294104140, 0.60725293651701, 0.60725293538591, 0.60725293510314,
0.60725293503245, 0.60725293501477, 0.60725293501035, 0.60725293500925,
0.60725293500897, 0.60725293500890, 0.60725293500889, 0.60725293500888
]
radians = lambda degrees: (degrees * pi) / 180
 
def Cordic(alpha, n, c_cos, c_sin):
i, ix, sigma = 0, 0, 0
kn, x, y, atn, t, theta = 0.0, 0.0, 0.0, 0.0, 0.0, 0.0
pow2 = 1.0
newsgn = 1 if floor(alpha / (2.0 * pi)) % 2 == 1 else -1
if alpha < -pi/2.0 or alpha > pi/2.0:
if alpha < 0:
Cordic(alpha + pi, n, x, y)
else:
Cordic(alpha - pi, n, x, y)
c_cos = x * newsgn
c_sin = y * newsgn
return c_cos, c_sin
ix = n - 1
if ix > 23: ix = 23
kn = kvalues[ix]
x = 1
y = 0
for i in range(n):
atn = angles[i]
sigma = 1 if theta < alpha else -1
theta += sigma * atn
t = x
x -= sigma * y * pow2
y += sigma * t * pow2
pow2 /= 2.0
c_cos = x * kn
c_sin = y * kn
return c_cos, c_sin
 
def main():
i, th = 0, 0
thr, c_cos, c_sin = 0.0, 0.0, 0.0
angles = [-9.0, 0.0, 1.5, 6.0]
print(" x sin(x) diff. sine cos(x) diff. cosine")
for th in range(-90,91,15):
thr = radians(th)
c_cos, c_sin = Cordic(thr, 24, c_cos, c_sin)
sin_diff = c_sin - sin(thr)
cos_diff = c_cos - cos(thr)
print(f"{th:+03}.0° {c_sin:+.8f} ({sin_diff:+.8f}) {c_cos:+.8f} ({cos_diff:+.8f})")
print("\nx(rads) sin(x) diff. sine cos(x) diff. cosine")
for i in range(4):
thr = angles[i]
c_cos, c_sin = Cordic(thr, 24, c_cos, c_sin)
sin_diff = c_sin - sin(thr)
cos_diff = c_cos - cos(thr)
print(f"{thr:+4.1f} {c_sin:+.8f} ({c_sin - sin(thr):+.8f}) {c_cos:+.8f} ({c_cos - cos(thr):+.8f})")
 
if __name__ in "__main__":
main()
</syntaxhighlight>
{{out}}
<pre>
x sin(x) diff. sine cos(x) diff. cosine
-90.0° -1.00000000 (+0.00000000) -0.00000007 (-0.00000007)
-75.0° -0.96592585 (-0.00000003) +0.25881895 (-0.00000009)
-60.0° -0.86602545 (-0.00000005) +0.49999992 (-0.00000008)
-45.0° -0.70710684 (-0.00000006) +0.70710672 (-0.00000006)
-30.0° -0.49999992 (+0.00000008) +0.86602545 (+0.00000005)
-15.0° -0.25881895 (+0.00000009) +0.96592585 (+0.00000003)
+00.0° +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+15.0° +0.25881895 (-0.00000009) +0.96592585 (+0.00000003)
+30.0° +0.49999992 (-0.00000008) +0.86602545 (+0.00000005)
+45.0° +0.70710684 (+0.00000006) +0.70710672 (-0.00000006)
+60.0° +0.86602545 (+0.00000005) +0.49999992 (-0.00000008)
+75.0° +0.96592585 (+0.00000003) +0.25881895 (-0.00000009)
+90.0° +1.00000000 (-0.00000000) -0.00000007 (-0.00000007)
 
x(rads) sin(x) diff. sine cos(x) diff. cosine
-9.0 -0.00000000 (+0.41211849) -0.00000000 (+0.91113026)
+0.0 +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+1.5 +0.99749499 (+0.00000000) +0.07073719 (-0.00000002)
+6.0 -0.00000000 (+0.27941550) -0.00000000 (-0.96017029)
</pre>
=={{header|Raku}}==
{{trans|XPL0}}
<syntaxhighlight lang="raku" line># 2023082 Raku programming solution
 
sub CORDIC ($A is copy) {
 
my (\Ten, $K, $X, $Y) = ( 1, * * 1/10 ... * )[^6], 0, 1, 0;
 
my \Tbl = < 7.853981633974480e-1 9.966865249116200e-2 9.999666686665240e-3
9.999996666668670e-4 9.999999966666670e-5 9.999999999666670e-6 0.0>;
 
while $A > 1e-5 {
$K++ while $A < Tbl[$K];
$A -= Tbl[$K];
($X,$Y) = $X - Ten[$K]*$Y, $Y + Ten[$K]*$X;
}
return $X, sqrt($X*$X + $Y*$Y)
}
 
say "Angle CORDIC Cosine Error";
for <-9 0 1.5 6> {
my \result = [/] CORDIC .abs;
printf "% 2.1f "~"% 2.8f " x 3~"\n", $_, result, .cos, .cos - result
} </syntaxhighlight>
{{out}} Same as XPL0 or you may [https://ato.pxeger.com/run?1=XZLfTsIwFMYTL_cUX5aZ8GcrK4OyBSQh6IXhwsRgAgE0YIYSccNuixKCL-KNF_pSPo2nKwvGNm2633fO2Tmn_fiS86fs8_M7S5eO_3Nyk2QL9K-uzy_7KFk9rBLcx5ttGTvDAPC8RWk6DCMb1oDWiNa4jDOUwG1UaPIad8EYo2N5citmNlxbaW678J8OF2vy6KDF_KYX-Fx4XtBqNHw3dDgCFgjhi2a9EXAu6i7BuoIBYSUIJRH0oML9HbmRNlOGLTJqFLDACjaPUGMFBVzmdnWOr4-rdQiqvQuuzHfFr6xBtXpUO6BCJtZg1i7kHpyz_7BEPdItskZwQK1TasUaq86hegSj3GOvNhmmmYzy7iYvMqUQJJOtNSa_srE3jGS-hdmLHigTGofr0qMfJ6soPHxcSBlLs20sY4mOE8AFZ02I7qEmdRsyTLJ1SvlNarMiEpsvkjydjVxF6RLmKeqMLwHzPT_66og3eO_mNDKpkjsbOo4Ndh8neqdyNTT2-nUdHlnx2H4B Attempt This Online!]
 
=={{header|RPL}}==
Line 577 ⟶ 1,192:
2: "6 constants in memory."
1: { -.91113026188 1 7.07372016661E-2 .960170286655 }
</pre>
 
=={{header|Rust}}==
{{trans|C}}
<syntaxhighlight lang="Rust">
use std::f64::consts::PI;
 
// Constants for angles and kvalues
const ANGLES: [f64; 28] = [
0.78539816339745, 0.46364760900081, 0.24497866312686, 0.12435499454676,
0.06241880999596, 0.03123983343027, 0.01562372862048, 0.00781234106010,
0.00390623013197, 0.00195312251648, 0.00097656218956, 0.00048828121119,
0.00024414062015, 0.00012207031189, 0.00006103515617, 0.00003051757812,
0.00001525878906, 0.00000762939453, 0.00000381469727, 0.00000190734863,
0.00000095367432, 0.00000047683716, 0.00000023841858, 0.00000011920929,
0.00000005960464, 0.00000002980232, 0.00000001490116, 0.00000000745058,
];
 
const KVALUES: [f64; 24] = [
0.70710678118655, 0.63245553203368, 0.61357199107790, 0.60883391251775,
0.60764825625617, 0.60735177014130, 0.60727764409353, 0.60725911229889,
0.60725447933256, 0.60725332108988, 0.60725303152913, 0.60725295913894,
0.60725294104140, 0.60725293651701, 0.60725293538591, 0.60725293510314,
0.60725293503245, 0.60725293501477, 0.60725293501035, 0.60725293500925,
0.60725293500897, 0.60725293500890, 0.60725293500889, 0.60725293500888,
];
 
// Function to convert degrees to radians
fn radians(degrees: f64) -> f64 {
degrees * PI / 180.0
}
 
// Cordic algorithm implementation
fn cordic(alpha: f64, n: usize, c_cos: &mut f64, c_sin: &mut f64) {
let mut theta = 0.0;
let mut pow2 = 1.0;
let mut x = 1.0;
let mut y = 0.0;
let newsgn = if (alpha / (2.0 * PI)).floor() as i32 % 2 == 1 { 1.0 } else { -1.0 };
 
if alpha < -PI / 2.0 || alpha > PI / 2.0 {
if alpha < 0.0 {
cordic(alpha + PI, n, &mut x, &mut y);
} else {
cordic(alpha - PI, n, &mut x, &mut y);
}
*c_cos = x * newsgn;
*c_sin = y * newsgn;
return;
}
 
let ix = if n - 1 > 23 { 23 } else { n - 1 };
let kn = KVALUES[ix];
 
for i in 0..n {
let atn = ANGLES[i];
let sigma = if theta < alpha { 1.0 } else { -1.0 };
theta += sigma * atn;
let t = x;
x -= sigma * y * pow2;
y += sigma * t * pow2;
pow2 /= 2.0;
}
 
*c_cos = x * kn;
*c_sin = y * kn;
}
 
fn main() {
let mut c_cos: f64;
let mut c_sin: f64;
let test_angles = [-9.0, 0.0, 1.5, 6.0];
 
println!(" x sin(x) diff. sine cos(x) diff. cosine");
for th in (-90..=90).step_by(15) {
let thr = radians(th as f64);
c_cos = 0.0;
c_sin = 0.0;
cordic(thr, 24, &mut c_cos, &mut c_sin);
println!("{:+03}.0° {:+.8} ({:+.8}) {:+.8} ({:+.8})", th, c_sin, c_sin - thr.sin(), c_cos, c_cos - thr.cos());
}
 
println!("\nx(rads) sin(x) diff. sine cos(x) diff. cosine");
for &angle in &test_angles {
let thr = angle;
c_cos = 0.0;
c_sin = 0.0;
cordic(thr, 24, &mut c_cos, &mut c_sin);
println!("{:+4.1} {:+.8} ({:+.8}) {:+.8} ({:+.8})", thr, c_sin, c_sin - thr.sin(), c_cos, c_cos - thr.cos());
}
}
</syntaxhighlight>
{{out}}
<pre>
x sin(x) diff. sine cos(x) diff. cosine
-90.0° -1.00000000 (+0.00000000) -0.00000007 (-0.00000007)
-75.0° -0.96592585 (-0.00000003) +0.25881895 (-0.00000009)
-60.0° -0.86602545 (-0.00000005) +0.49999992 (-0.00000008)
-45.0° -0.70710684 (-0.00000006) +0.70710672 (-0.00000006)
-30.0° -0.49999992 (+0.00000008) +0.86602545 (+0.00000005)
-15.0° -0.25881895 (+0.00000009) +0.96592585 (+0.00000003)
+00.0° +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+15.0° +0.25881895 (-0.00000009) +0.96592585 (+0.00000003)
+30.0° +0.49999992 (-0.00000008) +0.86602545 (+0.00000005)
+45.0° +0.70710684 (+0.00000006) +0.70710672 (-0.00000006)
+60.0° +0.86602545 (+0.00000005) +0.49999992 (-0.00000008)
+75.0° +0.96592585 (+0.00000003) +0.25881895 (-0.00000009)
+90.0° +1.00000000 (-0.00000000) -0.00000007 (-0.00000007)
 
x(rads) sin(x) diff. sine cos(x) diff. cosine
-9.0 -0.41211842 (+0.00000006) -0.91113029 (-0.00000003)
+0.0 +0.00000007 (+0.00000007) +1.00000000 (-0.00000000)
+1.5 +0.99749499 (+0.00000000) +0.07073719 (-0.00000002)
+6.0 -0.27941552 (-0.00000002) +0.96017028 (-0.00000001)
 
</pre>
 
=={{header|Sidef}}==
{{trans|Perl}}
<syntaxhighlight lang="ruby">func CORDIC(a) {
var (k, x, y) = (0, 1, 0)
var PoT = %n(1, 0.1, 0.01, 0.001, 0.0001, 0.00001)
var Tbl = %n(7.853981633e-1 9.966865249e-2 9.999666686e-3 9.999996666e-4 9.999999966e-5 9.999999999e-6 0)
 
for (true; a > 1e-5; a -= Tbl[k]) {
while (a < Tbl[k]) { ++k }
(x,y) = (x - PoT[k]*y, y + PoT[k]*x)
}
x / hypot(x,y)
}
 
print("Angle CORDIC Cosine Error\n")
for angle in (%n(-9 0 1.5 6)) {
var cordic = CORDIC(angle.abs)
printf("%4.1f %12.8f %12.8f %12.8f\n", angle, cordic, cos(angle), cos(angle) - cordic)
}</syntaxhighlight>
{{out}}
<pre>
Angle CORDIC Cosine Error
-9.0 -0.91112769 -0.91113026 -0.00000257
0.0 1.00000000 1.00000000 0.00000000
1.5 0.07073880 0.07073720 -0.00000160
6.0 0.96016761 0.96017029 0.00000268
</pre>
 
Line 582 ⟶ 1,340:
{{libheader|Wren-fmt}}
This is based on the Python example in the Wikipedia article except that the constants have been pre-calculated and the angles are adjusted, where necessary, to bring them within the [-π/2, π/2] range. The number of iterations is taken as 24 to try and match the Julia output which it does except for sin(0) which, curiously, has the same magnitude but a different sign.
<syntaxhighlight lang="ecmascriptwren">import "./fmt" for Fmt
 
// The following are pre-computed to avoid using atan and sqrt functions.
7,795

edits