Commit 4f3f4226 authored by Kevin Bracey's avatar Kevin Bracey
Browse files

* Added C99 functions floorf, ceilf, [l]rint[f], trunc[f], [l]round[f],...

* Added C99 functions floorf, ceilf, [l]rint[f], trunc[f], [l]round[f], nearbyint[f] and remainder[f].

* Added macros (only) for sinf, cosf, atanf.
* Requires cc 5.43

Version 5.39. Tagged as 'RISC_OSLib-5_39'
parent 84335a9f
......@@ -12,14 +12,14 @@
GBLS Module_HelpVersion
GBLS Module_ComponentName
GBLS Module_ComponentPath
Module_MajorVersion SETS "5.38"
Module_Version SETA 538
Module_MajorVersion SETS "5.39"
Module_Version SETA 539
Module_MinorVersion SETS ""
Module_Date SETS "26 Mar 2002"
Module_ApplicationDate2 SETS "26-Mar-02"
Module_ApplicationDate4 SETS "26-Mar-2002"
Module_Date SETS "11 Apr 2002"
Module_ApplicationDate2 SETS "11-Apr-02"
Module_ApplicationDate4 SETS "11-Apr-2002"
Module_ComponentName SETS "RISC_OSLib"
Module_ComponentPath SETS "RiscOS/Sources/Lib/RISC_OSLib"
Module_FullVersion SETS "5.38"
Module_HelpVersion SETS "5.38 (26 Mar 2002)"
Module_FullVersion SETS "5.39"
Module_HelpVersion SETS "5.39 (11 Apr 2002)"
END
/* (5.38)
/* (5.39)
*
* This file is automatically maintained by srccommit, do not edit manually.
* Last processed by srccommit version: 1.62.
*
*/
#define Module_MajorVersion_CMHG 5.38
#define Module_MajorVersion_CMHG 5.39
#define Module_MinorVersion_CMHG
#define Module_Date_CMHG 26 Mar 2002
#define Module_Date_CMHG 11 Apr 2002
#define Module_MajorVersion "5.38"
#define Module_Version 538
#define Module_MajorVersion "5.39"
#define Module_Version 539
#define Module_MinorVersion ""
#define Module_Date "26 Mar 2002"
#define Module_Date "11 Apr 2002"
#define Module_ApplicationDate2 "26-Mar-02"
#define Module_ApplicationDate4 "26-Mar-2002"
#define Module_ApplicationDate2 "11-Apr-02"
#define Module_ApplicationDate4 "11-Apr-2002"
#define Module_ComponentName "RISC_OSLib"
#define Module_ComponentPath "RiscOS/Sources/Lib/RISC_OSLib"
#define Module_FullVersion "5.38"
#define Module_HelpVersion "5.38 (26 Mar 2002)"
#define Module_LibraryVersionInfo "5:38"
#define Module_FullVersion "5.39"
#define Module_HelpVersion "5.39 (11 Apr 2002)"
#define Module_LibraryVersionInfo "5:39"
......@@ -45,6 +45,7 @@
#ifndef NO_FLOATING_POINT
#include <math.h> /* for forward references */
#include <fenv.h>
/* On the ARM, this has moved into the library's static data area */
/* so that it still works with the Shared C Library module. */
......@@ -615,79 +616,85 @@ double atan2(double y, double x)
return y;
}
/* Use the Moler-Morrison algorithm */
#ifndef HOST_HAS_TRIG
double hypot(double x, double y)
{
double p, q, r, s;
#pragma STDC FENV_ACCESS ON
int excepts;
double p, q, r;
fenv_t env;
p = fabs(x);
q = fabs(y);
/* Ensure p <= q (putting NaNs in p and Infs in q). We avoid calling
* functions like isnan() and isinf(), with their unwanted niceties
* of not triggering signalling NaNs. If we can avoid calling any
* functions, we save a lot of function call overhead.
*
* This ends up compiling neatly to:
*
* CMF p, q
* BGT swap
* CMFVS p, p
* BVC cont
* swap MVFD r, p
* MVFD p, q
* MVFD q, r
* cont
*
*/
if (isgreater(p,q) || isunordered(p,q) && isunordered(p,p))
if (isunordered(p, q))
{
if (p == INFINITY || q == INFINITY)
return INFINITY;
else
return p + q;
}
if (isless(p, q))
{
r = p; p = q; q = r;
}
/* More cunning ordering to only do two comparisons */
if (q == INFINITY || p == 0.0) return q;
if (isunordered(p, 0.0)) return p;
if (p == 0 || q == INFINITY)
return q;
for (;;)
feholdexcept(&env);
r = q / p;
p = p * sqrt(1 + r * r);
excepts = fetestexcept(FE_OVERFLOW);
fesetenv(&env);
if (excepts)
{
r = q / p;
r = r * r;
s = r + 4.0;
if (!islessgreater(s, 4.0)) return p;
r = r / s;
p = p + 2.0 * r * p;
q = q * r;
errno = ERANGE;
return HUGE_VAL;
}
return p;
}
/* All over again, in floats */
float hypotf(float x, float y)
{
float p, q, r, s;
#pragma STDC FENV_ACCESS ON
int excepts;
float p, q; double r;
fenv_t env;
float sqrtf(float);
p = fabsf(x);
q = fabsf(y);
if (isgreater(p,q) || isunordered(p,q) && isunordered(p,p))
if (isunordered(p, q))
{
if (p == INFINITY || q == INFINITY)
return INFINITY;
else
return p + q;
}
if (isless(p, q))
{
r = p; p = q; q = r;
}
if (q == INFINITY || p == 0.0F) return q;
if (isunordered(p, 0.0F)) return p;
if (p == 0 || q == INFINITY)
return q;
for (;;)
feholdexcept(&env);
r = (double) q / p;
p = (float) (p * sqrt(1 + r * r));
excepts = fetestexcept(FE_OVERFLOW);
fesetenv(&env);
if (excepts)
{
r = q / p;
r = r * r;
s = r + 4.0F;
if (!islessgreater(s, 4.0F)) return p;
r = r / s;
p = p + 2.0F * r * p;
q = q * r;
errno = ERANGE;
return HUGE_VAL;
}
return p;
}
#endif
/* Use the macros, which expand to inline ABS instructions */
double (fabs)(double x) { return fabs(x); }
......@@ -819,7 +826,7 @@ double fmod(double x, double y)
return r;
}
double floor(double d)
double (floor)(double d)
{
/* round x down to an integer towards minus infinity. */
fp_number x;
......@@ -846,7 +853,7 @@ double floor(double d)
else return x.d;
}
double ceil(double d)
double (ceil)(double d)
{
/* round x up to an integer towards plus infinity. */
fp_number x;
......@@ -873,6 +880,24 @@ double ceil(double d)
else return x.d;
}
/*double (ceil)(double x) { return ceil(x); }*/
float (ceilf)(float x) { return ceilf(x); }
/*double (floor)(double x) { return floor(x); }*/
float (floorf)(float x) { return floorf(x); }
double (rint)(double x) { return rint(x); }
float (rintf)(float x) { return rintf(x); }
long int (lrint)(double x) { return lrint(x); }
long int (lrintf)(float x) { return lrintf(x); }
double (trunc)(double x) { return trunc(x); }
float (truncf)(float x) { return truncf(x); }
double (atan)(double x) { return atan(x); }
/*float (atanf)(float x) { return atanf(x); }*/
double (sin)(double x) { return sin(x); }
/*float (sinf)(float x) { return sinf(x); }*/
double (cos)(double x) { return cos(x); }
/*float (cosf)(float x) { return cosf(x); }*/
double modf(double value, double *iptr)
{
/* splits value into integral part & fraction (both same sign) */
......
......@@ -15,10 +15,10 @@
#pragma force_top_level
#pragma include_only_once
/* math.h: ANSI 'C' (X3J11 Oct 88) library header, section 4.5 */
/* math.h: ISO 'C' (9899:1999) library header, section 7.22 */
/* Copyright (C) Codemist Ltd. */
/* Copyright (C) Acorn Computers Ltd. 1991 */
/* version 0.03 */
/* version 0.04 */
#ifndef __math_h
#define __math_h
......@@ -120,7 +120,9 @@ extern double atan2(double /*x*/, double /*y*/);
/* Returns: the arc tangent of y/x, in the range -Pi to Pi. */
extern double __d_atan(double);
extern __caller_narrow float __r_atan(float);
#define atan(x) __d_atan(x)
#define atanf(x) __r_atan(x)
extern double cos(double /*x*/);
/* computes the cosine of x (measured in radians). A large magnitude */
......@@ -132,9 +134,13 @@ extern double sin(double /*x*/);
/* Returns: the sine value. */
extern double __d_sin(double);
extern __caller_narrow float __r_sin(float);
extern double __d_cos(double);
extern __caller_narrow float __r_cos(float);
#define sin(x) __d_sin(x)
#define sinf(x) __r_sin(x)
#define cos(x) __d_cos(x)
#define cosf(x) __r_cos(x)
extern double tan(double /*x*/);
/* computes the tangent of x (measured in radians). A large magnitude */
......@@ -208,27 +214,83 @@ extern double sqrt(double /*x*/);
/* if the argument is negative, and -HUGE_VAL returned. */
/* Returns: the value of the square root. */
extern double ceil(double /*x*/);
/* computes the smallest integer not less than x. */
/* Returns: the smallest integer not less than x, expressed as a double. */
extern double fabs(double /*x*/);
extern float fabsf(float /*x*/);
/* computes the absolute value of the floating-point number x. */
/* Returns: the absolute value of x. */
extern double floor(double /*d*/);
extern double ceil(double /*x*/);
extern float ceilf(float /*x*/);
/* computes the smallest integer not less than x. */
/* Returns: the smallest integer not less than x, expressed as a double. */
extern double floor(double /*x*/);
extern float floorf(float /*x*/);
/* computes the largest integer not greater than x. */
/* Returns: the largest integer not greater than x, expressed as a double */
extern double nearbyint(double /*x*/);
extern float nearbyintf(float /*x*/);
/* rounds its argument to an integer value, using the current rounding */
/* direction. Does not raise the inexact exception. */
/* Returns: the rounded integer value. */
extern double rint(double /*x*/);
extern float rintf(float /*x*/);
/* rounds its argument to an integer value, using the current rounding */
/* direction. Raises "inexact" if the result differs from the argument. */
/* Returns: the rounded integer value. */
extern long int lrint(double /*x*/);
extern long int lrintf(float /*x*/);
/* rounds its argument to an integer value, using the current rounding */
/* direction. Raises "inexact" if the result differs from the argument. */
/* Returns: the rounded integer value. */
extern double round(double /*x*/);
extern float roundf(float /*x*/);
/* rounds its argument to the nearest integer value, rounding halfway */
/* cases away from zero. */
/* Returns: the rounded integer value. */
extern long int lround(double /*x*/);
extern long int lroundf(float /*x*/);
/* rounds its argument to the nearest integer value, rounding halfway */
/* cases away from zero. */
/* Returns: the rounded integer value. */
extern double trunc(double /*x*/);
extern float truncf(float /*x*/);
/* rounds its argument to the integer value, nearest to but no larger in */
/* magnitude than the argument. */
/* Returns: the truncated integer value. */
extern double fmod(double /*x*/, double /*y*/);
/* computes the floating-point remainder of x/y. */
/* Returns: the value x - i * y, for some integer i such that, if y is */
/* nonzero, the result has the same sign as x and magnitude */
/* less than the magnitude of y. If y is zero, a domain error */
/* occurs and -HUGE_VAL is returned. */
extern double remainder(double /*x*/, double /*y*/);
extern float remainderf(float /*x*/, float /*y*/);
/* computes the remainder x REM y required by IEEE 754 */
/* Returns: x REM y */
extern double __d_abs(double);
extern float __r_abs(float);
extern double __d_floor(double);
extern double __d_ceil(double);
extern double __d_trunc(double);
extern double __d_rint(double);
extern long int __d_lrint(double);
extern float __r_abs(float);
extern __caller_narrow float __r_floor(float);
extern __caller_narrow float __r_ceil(float);
extern __caller_narrow float __r_trunc(float);
extern __caller_narrow float __r_rint(float);
extern __caller_narrow long int __r_lrint(float);
#define fabs(x) __d_abs(x)
#define fabsf(x) __r_abs(x)
#define floor(x) __d_floor(x)
#define floorf(x) __r_floor(x)
#define ceil(x) __d_ceil(x)
#define ceilf(x) __r_ceil(x)
#define trunc(x) __d_trunc(x)
#define truncf(x) __r_trunc(x)
#define rint(x) __d_rint(x)
#define rintf(x) __r_rint(x)
#define lrint(x) __d_lrint(x)
#define lrintf(x) __r_lrint(x)
#pragma force_fpargs_in_regs
extern double copysign(double /*x*/, double /*y*/);
......
......@@ -58,10 +58,13 @@ $Label MOV r1, #0
MEND
MACRO
$Label ReEnableFPInterrupts
$Label ReEnableFPInterrupts $nocheck
; Reinstates the exception mask state which prevailed before the call
; to DisableFPInterrupts; sets r1 to the current fp flags.
$Label RFS r1
$Label
[ "$nocheck"=""
RFS r1
]
WFS ip
MEND
......@@ -163,17 +166,25 @@ $Label RFS r1
EXPORT fmax
EXPORT fmaxf
EXPORT sin
EXPORT cos
EXPORT nearbyint
EXPORT nearbyintf
EXPORT round
EXPORT roundf
EXPORT lround
EXPORT lroundf
EXPORT remainder
EXPORT remainderf
EXPORT exp
EXPORT log10
EXPORT log
EXPORT sqrt
EXPORT sqrtf
EXPORT tan
EXPORT atan
EXPORT asin
EXPORT acos
EXPORT pow
EXPORT hypot
EXPORT hypotf
GET clib.s.h_signal
......@@ -1145,43 +1156,110 @@ fmaf
Return ,LinkNotStacked
]
sin
[ :LNOT: FloatingPointArgsInRegs
STMFD sp!, {r0, r1}
nearbyint
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0,r1}
LDFD f0, [sp], #8
]
SIND f0, f0
DisableFPInterrupts
RNDD f0, f0
ReEnableFPInterrupts nocheck
Return ,LinkNotStacked
[ {FALSE}
sinf
[ :LNOT: FloatingPointArgsInRegs
STMFD sp!, {r0, r1}
nearbyintf
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0,r1}
LDFD f0, [sp], #8
]
MVFS f0, f0
SINS f0, f0
DisableFPInterrupts
RNDS f0, f0
ReEnableFPInterrupts nocheck
Return ,LinkNotStacked
]
cos
[ :LNOT: FloatingPointArgsInRegs
STMFD sp!, {r0, r1}
round
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0,r1}
LDFD f0, [sp], #8
]
COSD f0, f0
CMF f0, #0
Return ,LinkNotStacked,EQ ; return +/-0 intact
DisableFPInterrupts ; just to prevent "inexact"
ABSMID f0, f0
ADFD f0, f0, #0.5
RNDDZ f0, f0
MNFMID f0, f0
ReEnableFPInterrupts nocheck
Return ,LinkNotStacked
[ {FALSE}
cosf
[ :LNOT: FloatingPointArgsInRegs
STMFD sp!, {r0, r1}
roundf
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0,r1}
LDFD f0, [sp], #8
]
MVFS f0, f0
COSS f0, f0
CMF f0, #0
Return ,LinkNotStacked,EQ ; return +/-0 intact
DisableFPInterrupts ; just to prevent "inexact"
ABSMIS f0, f0
ADFS f0, f0, #0.5
RNDSZ f0, f0
MNFMIS f0, f0
ReEnableFPInterrupts nocheck
Return ,LinkNotStacked
lround
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0,r1}
LDFD f0, [sp], #8
]
lround2 RFS ip
AND r2, ip, #&10 ; remember cumulative flag
BIC r1, ip, #&10 ; clear cumulative inexact flag
WFS r1
FIX r0, f0 ; may trap - that's okay
RFS r1
TST r1, #&10 ; was inexact?
BEQ %FT10
; If we get here, we know inexact traps are disabled, and result is inexact
CMF f0, #0
ADFGTD f0, f0, #0.5 ; no new exceptions possible
SUFMID f0, f0, #0.5
FIX r0, f0
FIXM r0, f0
FIXP r0, f0
FIXZ r0, f0
10 ORR r1, r1, r2 ; re-cumulate inexact flag
WFS r1
Return ,LinkNotStacked
lroundf
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0,r1}
LDFD f0, [sp], #8
]
MVFS f0, f0
B lround2
remainder
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0-r3}
LDFD f0, [sp], #8
LDFD f1, [sp], #8
]
RMFD f0, f0, f1
Return ,LinkNotStacked
remainderf
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0-r3}
LDFD f0, [sp], #8
LDFD f1, [sp], #8
]
MVFS f0, f0
MVFS f1, f1
RMFS f0, f0, f1
Return ,LinkNotStacked
]
exp
[ :LNOT: FloatingPointArgsInRegs
......@@ -1227,6 +1305,7 @@ huge_errorf
LDFS f0, huge_valf
B Set_errno
huge_valf
infinity
DCD &7F800000
negative_error
......@@ -1278,7 +1357,6 @@ sqrt
SQTD f0, f0
Return ,LinkNotStacked
[ {FALSE}
sqrtf
[ :LNOT: FloatingPointArgsInRegs
STMFD sp!, {r0, r1}
......@@ -1289,7 +1367,6 @@ sqrtf
BMI negative_errorf
SQTS f0, f0
Return ,LinkNotStacked
]
tan
[ :LNOT: FloatingPointArgsInRegs
......@@ -1322,14 +1399,6 @@ tanf
B huge_errorf
]
atan
[ :LNOT: FloatingPointArgsInRegs
STMFD sp!, {r0, r1}
LDFD f0, [sp], #8
]
ATND f0, f0
Return ,LinkNotStacked
asin
[ :LNOT: FloatingPointArgsInRegs
STMFD sp!, {r0, r1}
......@@ -1399,6 +1468,80 @@ POWFirstArgZero
Return ,LinkNotStacked, GE
B negative_error
hypot
[ FloatingPointArgsInRegs
ABSD f0, f0
ABSD f1, f1
|
BIC r0, r0, #&80000000 ; faster than using ABS later
BIC r2, r2, #&80000000
STMFD sp!, {r0-r3}
LDFD f0, [sp], #8
LDFD f1, [sp], #8
]
LDFS f3, infinity
CMF f0, f1
BVS hypot_nan
MVFMID f2, f0
MVFMID f0, f1
MVFMID f1, f2
CMF f0, #0
CMFNE f1, f3
MVFEQD f0, f1
Return ,LinkNotStacked, EQ
DisableFPInterrupts
DVFE f2, f1, f0 ; 0 <= result <= 1; may as well use extra precision
MUFE f2, f2, f2
ADFE f2, f2, #1 ; 1 <= result <= 2
SQTE f2, f2 ; 1 <= result <= sqrt(2)
MUFD f0, f0, f2 ; this may overflow, but only slightly
ReEnableFPInterrupts
TST r1, #4
Return ,LinkNotStacked, EQ
B huge_error
hypot_nan
CMF f0, f3
CMFNE f1, f3
MVFEQD f0, f3
ADFNED f0, f0, f1
Return ,LinkNotStacked
hypotf
[ :LNOT:FloatingPointArgsInRegs
STMFD sp!, {r0-r3}
LDFD f0, [sp], #8
LDFD f1, [sp], #8
]
ABSS f0, f0
ABSS f1, f1
LDFS f3, infinity
CMF f0, f1
BVS hypot_nan
MVFMIS f2, f0
MVFMIS f0, f1
MVFMIS f1, f2
CMF f0, #0
CMFNE f1, f3
MVFEQS f0, f1
Return ,LinkNotStacked, EQ