OpenCores
URL https://opencores.org/ocsvn/forwardcom/forwardcom/trunk

Subversion Repositories forwardcom

Compare Revisions

  • This comparison shows the changes necessary to convert path 
    /forwardcom/libraries
    from Rev 109 to Rev 110
    ↔ Reverse comparison
  •

Rev 109 → Rev 110

/sincos.as
0,0 → 1,157
/********************************* sin.as ***********************************
* Author: Agner Fog
* date created: 2018-03-29
* Last modified: 2021-04-25
* Version: 1.11
* Project: ForwardCom library math.li
* Description: sin, cos, and tan functions. Calculate in radians, double precision
* The argument x can be a scalar or a vector
* The return value will be a vector with the same length
* C declaration: double sin(double x);
* C declaration: double cos(double x);
* C declaration: double tan(double x);
* C declaration: struct {double s; double c;} sincos(double x);
*
* This code is adapted from C++ vector class library www.github.com/vectorclass
* Copyright 2018-2021 GNU General Public License http://www.gnu.org/licenses
*****************************************************************************/
 
 
// define constants
% M_2_PI = 0.636619772367581343076 // 2./pi
% P0sin = -1.66666666666666307295E-1 // polynomial coefficients for sin
% P1sin = 8.33333333332211858878E-3
% P2sin = -1.98412698295895385996E-4
% P3sin = 2.75573136213857245213E-6
% P4sin = -2.50507477628578072866E-8
% P5sin = 1.58962301576546568060E-10
 
% P0cos = 4.16666666666665929218E-2 // polynomial coefficients for cos
% P1cos = -1.38888888888730564116E-3
% P2cos = 2.48015872888517045348E-5
% P3cos = -2.75573141792967388112E-7
% P4cos = 2.08757008419747316778E-9
% P5cos = -1.13585365213876817300E-11
 
% DP1 = 7.853981554508209228515625E-1 // modulo pi/2 for extended precision
% DP2 = 7.94662735614792836714E-9 // correction for extended precision modular artithmetic
% DP3 = 3.06161699786838294307E-17 // correction for extended precision modular artithmetic
 
 
code section execute align = 4
 
public _sin: function, reguse = 0, 0x7BF
public _cos: function, reguse = 0, 0x7BF
public _sincos: function, reguse = 0, 0x7BF
public _tan: function, reguse = 0, 0x7BF
 
// common entry for sin and sincos functions
_sin function
_sincos:
 
/* registers:
v0 = x
v1 = abs(x)
v1 = quadrant
v10 = abs(x) reduced modulo pi/2
v2 = v10^2
v3 = v10^4
v4 = v10^8
v5 = v10^3
v6 = unused (vacant flag for calling function)
v7 = temp
v8 = sin
v9 = cos
*/
 
// Find quadrant:
// 0 - pi/4 => 0
// pi/4 - 3*pi/4 => 1
// 3*pi/4 - 5*pi/4 => 2
// 5*pi/4 - 7*pi/4 => 3
// 7*pi/4 - 8*pi/4 => 4
 
double v1 = clear_bit(v0, 63) // abs(x)
double v4 = v1 * M_2_PI
double v4 = round(v4, 0) // round to integer
 
// reduce modulo pi/2, with extended precision
// x = ((xa - y * DP1) - y * DP2) - y * DP3;
double v10 = v4 * (-DP1*2) + v1
double v10 = v4 * (-DP2*2) + v10
double v10 = v4 * (-DP3*2) + v10
double v5 = !(v4 > ((1 << 51) + 0.0)) // check for loss of precision and overflow, but not NAN
double v1 = v4 + ((1 << 52) + 0.0) // add magic number 2^52 to get integer into lowest bit
double v10 = v5 ? v10 : 0 // zero if out of range. result will be -1, 0, or 1
 
// Expansion of sin and cos, valid for -pi/4 <= x <= pi/4
double v2 = v10 * v10 // x^2
double v3 = v2 * v2 // x^4
double v4 = v3 * v3 // x^8
 
// calculate polynomial P5sin*x2^5 + P4sin*x2^4 + P3sin*x2^3 + P2sin*x2^2 + P1sin*x2 + P0sin
// = (p2+p3*x2)*x4 + ((p4+p5*x2)*x8 + (p0+p1*x2));
 
double v5 = replace(v10, P0sin) // broadcast to same length as x
double v5 = v2 * P1sin + v5
double v7 = replace(v10, P4sin)
double v7 = v2 * P5sin + v7
double v8 = replace(v10, P2sin)
double v8 = v2 * P3sin + v8
double v7 = v7 * v4 + v5
double v8 = v8 * v3 + v7
 
// calculate polynomial P5cos*x2^5 + P4cos*x2^4 + P3cos*x2^3 + P2cos*x2^2 + P1cos*x2 + P0cos
// = (p2+p3*x2)*x4 + ((p4+p5*x2)*x8 + (p0+p1*x2));
double v5 = replace(v10, P0cos)
double v5 = v2 * P1cos + v5
double v7 = replace(v10, P4cos)
double v7 = v2 * P5cos + v7
double v9 = replace(v10, P2cos)
double v9 = v2 * P3cos + v9
double v7 = v7 * v4 + v5
double v9 = v9 * v3 + v7
 
// s = x + (x * x2) * s;
double v5 = v10 * v2
double v8 = v8 * v5 + v10
// c = 1.0 - x2 * 0.5 + (x2 * x2) * c;
double v9 = v9 * v3
double v9 = v2 * (-0.5) + v9
double v9 = 1.0 + v9
 
// swap sin and cos if odd quadrant
double v3 = v1 ? v9 : v8 // sin
double v4 = v1 ? v8 : v9 // cos
 
// get sign of sin
int64 v5 = v1 << 62 // get bit 1 into sign bit, x modulo pi/2 = 2 or 3
int64 v5 ^= v0 // toggle with sign of original x
int64 v5 = and(v5, 1 << 63) // isolate sign bit
double v0 = v3 ^ v5 // apply sign bit to sin
 
// get sign of cos
int64 v1 = v1 + 1 // change sign when x modulo pi/2 = 1 or 2
int64 v1 = v1 << 62 // get bit 1 into sign bit
int64 v1 = and(v1, 1 << 63) // isolate sign bit
double v1 = v4 ^ v1 // apply sign bit to cos
 
// return sin in v0, cos in v1
return
_sin end
 
// cosine function
_cos function
call _sincos
double v0 = v1 // cos is in v1
return
_cos end
 
// tangent function
_tan function
call _sincos
double v0 = v0 / v1 // tan(x) = sin(x)/cos(x)
return
_tan end
 
code end

powered by: WebSVN 2.1.0

© copyright 1999-2024 OpenCores.org, equivalent to Oliscience, all rights reserved. OpenCores®, registered trademark.