Sunday, February 18, 2024

Casio fx-CG 50 CORDIC Simulation: Approximating Sine and Cosine of Angles

Casio fx-CG 50 CORDIC Simulation:  Approximating Sine and Cosine of Angles

Introduction - How Computers Calculate Mathematical Functions

First developed by Jack E. Volder, the Coordinate Rotation Digital Computer, better known as CORDIC, is an algorithm used to calculate many mathematical functions, including trigonometric functions, logarithms, exponentials, and hyperbolic functions.     CORDIC is a fundamental algorithm, which variants of CORDIC are used in computers and calculators.

Today's focus will be a calculating sines and cosines of angles.  The steps will be detailed in the next sections.

CORDIC1:  Series of Arctangents

Let Θ be the angle.   The first step is to build Θ as additions and subtractions of terms of arctan(1 ÷ (2^n)), starting at n = 0 and stopping at the required accuracy.  Let A be the approximation.

Examples:

Θ = 45° requires only one term:

A = 45° = arctan 1

Θ ≈ 71.56505118°

A = 71.56505118° = arctan 1 + arctan 1/2

Θ ≈ 32.47119229°

A = 32.47119229° = arctan 1 - arctan 1/2 + arctan 1/4

A series to recreate Θ = 30° takes 32 terms of ±arctan (1 ÷ 2^n) from n = 0 to n = 31.   (8 decimal places)  This is known as coordinate rotation.

The program CORDIC1 determines the number of terms need to created to obtain Θ.   Accuracy in this code is set to 5 decimal places, but can be adjusted (see the line in blue).

Casio fx-CG 50 Program Code:  CORDIC1

Deg

"ANGLE"?->Θ

0->I

0->A

Lbl 0

tan^-1 (2^(-I))->T

If A<Θ

Then

A+T->A

Else

A-T->A

IfEnd

I+1->I

Abs (A-Θ)>1×10^-5=>Goto 0

ClrText

Blue Locate 1,3,"Θ: "

Blue Locate 5,3,Θ

Red Locate 1,4,"A: "

Red Locate 5,4,A

Green Locate 1,7,I

Notes:

-> is the store arrow →

=> is the jump command ⇒

For reference:

arctan 1 = 45°

arctan 1/2 ≈ 26.56505118°

arctan 1/4 ≈ 14.03624347°

arctan 1/8 ≈ 7.125016349°

arctan 1/16 ≈ 3.576334375°

arctan 1/32 ≈ 1.789910608°

CORDIC2:  Calculating Sine and Cosine

Imagine the coordinate (cos Θ, sin Θ) on a unit circle.  A unit circle is a circle with radius of length 1 and center located at the origin (0,0).

Set the initial angle at 0°.  Then x = cos 0° = 1 and y = sin 0° = 0.   The initial vector is set to be [ [ cos 0° ] [ sin 0° ] ] = [ [ 1 ] [ 0 ] ].

Approximate Θ in terms of arctangent (1 ÷ (2^n)), starting at n = 0.

Direction of Rotation

Set σ as the direction of rotation.

A rotation is positive if we add arctan(1 ÷ (2^n)) to A.    For a positive rotation, set σ_i = +1.

A rotation is negative if we subtract arctan(1 ÷ (2^n)) to A.    For a negative rotation, set σ_i = -1.

For example:

Θ ≈ 32.47119229°

A = 32.47119229° = arctan 1 - arctan 1/2 + arctan 1/4

Then:  σ_0 = 1 (positive rotation), σ_1 = -1 (negative rotation), σ_2 =  1 (positive rotation).

Multiplying Factor

The number of iterations is also used to determine the required multiplication factor:

n = Π( 1 ÷ √(1 + 2^(-2 × I)), I = 0 to I = terms needed - 1)

For the example:

Θ ≈ 32.47119229°

Needed 3 iterations, i_0 to i_2.

Then:

n = 1 ÷ √(1 + 2^(-2 × 0)) × 1 ÷ √(1 + 2^(-2 × 1)) × 1 ÷ √(1 + 2^(-2 × 2))

= 4 × √170 ÷ 85

≈ 0.6135719911

Calculating the Next Iteration

The next iteration for x and y are:

x_i+1 = x_i - 2^(-i) × σ_i × y_i

y_i+1 = 2^(-i) × σ_i × x_i + y_i

When A is sufficiently near or equal to Θ, the cosine and sine are approximated as:

cos(A) ≈ n × x_final

sin(A) ≈ n × y_final

For the example:

Θ ≈ 32.47119229°

The approximated angle, in this case, happens to be exact angle, A = Θ.  And,

cos(A) ≈ 0.8436614877

sin(A) ≈ 0.5368754922

For more details, please check out resources in the Sources section.   I particularly like Oxford's A Very Short Introduction series.

Casio fx-CG 50 Program Code:  CORDIC2

This algorithm is adopted from the Python example (see Wikipedia article).

Deg

"ANGLE"?->Θ

0->I

0->A

1->X

0->Y

1->N

Lbl 0

tan^-1 (2^(-I))->T

If A<Θ

Then

A+T->A

1->S

Else

A-T->A

-1->S

IfEnd

N÷√(1+2^(-2*I))->N

X-S*2^(-I)*Y->P

Y+X*S*2^(-I)->Q

P->X

Q->Y

I+1->I

Abs (A-Θ)>1×10^-5=>Goto 0

X*N->X

Y*N->Y

ClrText

Blue Locate 1,3,"Θ: "

Blue Locate 5,3,Θ

Red Locate 1,4,"A: "

Red Locate 5,4,A

Black Locate 1,5,"cos :"

Black Locate 7,5,X

Black Locate 1,6,"sin :"

Black Locate 7,6,Y

Green Locate 1,7,I

CORDIC3:  Doing it without a Arctangent function

In reality, most of the time the arctangent function also has to be approximated.   There are many ways to approximate, which varying accuracy.   I used the fx-CG50's statistics mode to come up with a regression equation with the following lists:

x_list = sequence of 1÷(2^i) from i = 0 to i = 39

y_list = sequence of arctan(1÷(2^i)) from i = 0 to i = 39

Of the regression models the fx-CG50 offers, the best regression model is quartic regression (4th-order polynomial):

y ≈ 9.43597784 × x^4 - 22.116232 × x^3 + 0.40116676 × x^2 + 57.2790727 × x + 1.4558 × 10^-5

Remember that I am working with degree angle measurement.

Casio fx-CG 50 Program Code:  CORDIC3

Deg

"ANGLE"?->Θ

0->I

0->A

1->X

0->Y

1->N

Lbl 0

2^(-I)->K

9.43597784917955K^(4)-22.1162323028173K^(3)+

0.401166766193516K^2+57.2790727477367K+

1.45584715586876×10^-5->T

If A<Θ

Then

A+T->A

1->S

Else

A-T->A

-1->S

IfEnd

N÷√(1+2^(-2*I))->N

X-S*2^(-I)*Y->P

Y+X*S*2^(-I)->Q

P->X

Q->Y

I+1->I

Abs (A-Θ)>1×10^-5=>Goto 0

X*N->X

Y*N->Y

ClrText

Blue Locate 1,3,"Θ: "

Blue Locate 5,3,Θ

Red Locate 1,4,"A: "

Red Locate 5,4,A

Black Locate 1,5,"cos :"

Black Locate 7,5,X

Black Locate 1,6,"sin :"

Black Locate 7,6,Y

Green Locate 1,7,I

For the example:

Θ ≈ 32.47119229°

Approximating Θ within 8 decimal places (10^-5) yields these results:

A:  32.47118598

cos A:  0.8436683456

sin A:  0.5368647154

24 terms used

Sources

"CORDIC"  Wikipedia.   Last Edited January 11, 2024.   Accessed January 12, 2024.  https://en.wikipedia.org/wiki/CORDIC

Brummelen, Glen Van.  Trigonometry:  A Very Short Introduction  Oxford University Press: Oxford, United Kingdom.  2020.  ISBN 978-0-19-881431-3

Eddie

All original content copyright, © 2011-2024.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author.

HP Prime and Numworks: Plotting a Parametric Line of Motion

HP Prime and Numworks: Plotting a Parametric Line of Motion Plotting the Position of Motion This program draws a 2D motion plo...