mirror-games/glest-source
1
0
Fork 0
mirror of https://github.com/glest/glest-source.git synced 2025-08-12 03:14:00 +02:00
Code Issues Projects Releases Wiki Activity

Added streflop (standalone reproducible floating point library) layer to mega-glest (initial checkin only) and changed a few areas to use the library in linux

Browse Source
This commit is contained in:
Mark Vejvoda
2010-04-24 03:57:38 +00:00
parent a50b89e9cc
commit 2bfaa4d1d7
129 changed files with 15486 additions and 79 deletions
Download Patch File Download Diff File Expand all files Collapse all files

16
source/shared_lib/include/streflop/CMakeLists.txt Normal file
View File

@@ -0,0 +1,16 @@
AUX_SOURCE_DIRECTORY(libm/flt-32 libm_flt32_source)
SET(cxxflags "-w -O3 -I${CMAKE_CURRENT_SOURCE_DIR}/libm/headers")
if (NOT $ENV{CXX} MATCHES "icpc")
SET(cxxflags "${cxxflags} -mfpmath=sse -msse")
endif (NOT $ENV{CXX} MATCHES "icpc")
SET_SOURCE_FILES_PROPERTIES(${libm_flt32_source} PROPERTIES COMPILE_FLAGS "-DLIBM_COMPILING_FLT32 ${cxxflags}")
ADD_LIBRARY(streflop STATIC EXCLUDE_FROM_ALL
SMath.cpp
Random.cpp
streflopC.cpp
${libm_flt32_source}
)
set_target_properties(streflop PROPERTIES COMPILE_FLAGS "${PIC_FLAG}")
#TODO do not use -fPIC for streflop (decreases performance)

533
source/shared_lib/include/streflop/FPUSettings.h Normal file
View File

@@ -0,0 +1,533 @@
/*
streflop: STandalone REproducible FLOating-Point
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
Heavily relies on GNU Libm, itself depending on netlib fplibm, GNU MP, and IBM MP lib.
Uses SoftFloat too.
Please read the history and copyright information in the documentation provided with the source code
*/
/*
For reference, the layout of the MXCSR register:
FZ:RC:RC:PM:UM:OM:ZM:DM:IM:Rsvd:PE:UE:OE:ZE:DE:IE
15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 0
And the layout of the 387 FPU control word register:
Rsvd:Rsvd:Rsvd:X:RC:RC:PC:PC:Rsvd:Rsvd:PM:UM:OM:ZM:DM:IM
15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 0
Where:
Rsvd - Reserved
FZ - Flush to Zero
RC - Rounding Control
PM - Precision Mask
UM - Underflow Mask
OM - Overflow Mask
ZM - Zerodivide Mask
DM - Denormal Mask
IM - Invalid Mask
PE - Precision Exception
UE - Underflow Exception
OE - Overflow Exception
ZE - Zerodivide Exception
DE - Denormal Exception
IE - Invalid Exception
X - Infinity control (unused on 387 and higher)
PC - Precision Control
Source: Intel Architecture Software Development Manual, Volume 1, Basic Architecture
*/
// Included by the main streflop include file
// module broken apart for logical code separation
#ifndef STREFLOP_FPU_H
#define STREFLOP_FPU_H
// Can safely make the symbols from softfloat visible to user program, protected in namespace
#if defined(STREFLOP_SOFT)
#include "softfloat/softfloat.h"
#endif
namespace streflop {
// We do not use libm, so let's copy a few flags and C99 functions
// Give warning in case these flags would be defined already, this is indication
// of potential confusion!
#if defined(FE_INVALID) || defined(FE_DENORMAL) || defined(FE_DIVBYZERO) || defined(FE_OVERFLOW) || defined(FE_UNDERFLOW) || defined(FE_INEXACT) || defined(FE_DOWNWARD) || defined(FE_TONEAREST) || defined(FE_TOWARDZERO) || defined(FE_UPWARD)
#warning STREFLOP: FE_XXX flags were already defined and will be redefined! Check you do not use the system libm.
#undef FE_INVALID
#undef FE_DENORMAL
#undef FE_DIVBYZERO
#undef FE_OVERFLOW
#undef FE_UNDERFLOW
#undef FE_INEXACT
#undef FE_INEXACT
#undef FE_ALL_EXCEPT
#undef FE_DOWNWARD
#undef FE_TONEAREST
#undef FE_TOWARDZERO
#undef FE_UPWARD
#endif // defined(FE_INVALID) || ...
// Flags for FPU exceptions
enum FPU_Exceptions {
// Invalid operation. If not signaling, gives NaN instead
FE_INVALID = 0x0001,
#define FE_INVALID FE_INVALID
// Extension: for x86 and SSE
// Denormal operand. If not signaling, use denormal arithmetic as usual
FE_DENORMAL = 0x0002,
#define FE_DENORMAL FE_DENORMAL
// Division by zero. If not signaling, uses +/- infinity
FE_DIVBYZERO = 0x0004,
#define FE_DIVBYZERO FE_DIVBYZERO
// Overflow. If not signaling, round to nearest (including infinity) according to rounding mode
FE_OVERFLOW = 0x0008,
#define FE_OVERFLOW FE_OVERFLOW
// Underflow. If not signaling, use 0 instead
FE_UNDERFLOW = 0x0010,
#define FE_UNDERFLOW FE_UNDERFLOW
// Rounding was not exact (ex: sqrt(2) is never exact) or when overflow causes rounding
FE_INEXACT = 0x0020,
#define FE_INEXACT FE_INEXACT
// Combination of all the above
FE_ALL_EXCEPT = 0x003F
#define FE_ALL_EXCEPT FE_ALL_EXCEPT
};
// Flags for FPU rounding modes
enum FPU_RoundMode {
FE_TONEAREST = 0x0000,
#define FE_TONEAREST FE_TONEAREST
FE_DOWNWARD = 0x0400,
#define FE_DOWNWARD FE_DOWNWARD
FE_UPWARD = 0x0800,
#define FE_UPWARD FE_UPWARD
FE_TOWARDZERO = 0x0C00
#define FE_TOWARDZERO FE_TOWARDZERO
};
/* Note: SSE control word, bits 0..15
0->5: Run-time status flags
6: DAZ (denormals are zero, i.e. don't use denormals if bit is 1)
7->12: Exception flags, same meaning as for the x87 ones
13,14: Rounding flags, same meaning as for the x87 ones
15: Flush to zero (FTZ) for automatic handling of underflow (default is NO)
*/
// plan for portability
#if defined(_MSC_VER)
#define STREFLOP_FSTCW(cw) do { short tmp; __asm { fstcw tmp }; (cw) = tmp; } while (0)
#define STREFLOP_FLDCW(cw) do { short tmp = (cw); __asm { fclex }; __asm { fldcw tmp }; } while (0)
#define STREFLOP_STMXCSR(cw) do { int tmp; __asm { stmxcsr tmp }; (cw) = tmp; } while (0)
#define STREFLOP_LDMXCSR(cw) do { int tmp = (cw); __asm { ldmxcsr tmp }; } while (0)
#else // defined(_MSC_VER)
#define STREFLOP_FSTCW(cw) do { asm volatile ("fstcw %0" : "=m" (cw) : ); } while (0)
#define STREFLOP_FLDCW(cw) do { asm volatile ("fclex \n fldcw %0" : : "m" (cw)); } while (0)
#define STREFLOP_STMXCSR(cw) do { asm volatile ("stmxcsr %0" : "=m" (cw) : ); } while (0)
#define STREFLOP_LDMXCSR(cw) do { asm volatile ("ldmxcsr %0" : : "m" (cw) ); } while (0)
#endif // defined(_MSC_VER)
// Subset of all C99 functions
#if defined(STREFLOP_X87)
/// Raise exception for these flags
inline int feraiseexcept(FPU_Exceptions excepts) {
unsigned short fpu_mode;
STREFLOP_FSTCW(fpu_mode);
fpu_mode &= ~( excepts ); // generate error for selection
STREFLOP_FLDCW(fpu_mode);
return 0;
}
/// Clear exceptions for these flags
inline int feclearexcept(int excepts) {
unsigned short fpu_mode;
STREFLOP_FSTCW(fpu_mode);
fpu_mode |= excepts;
STREFLOP_FLDCW(fpu_mode);
return 0;
}
/// Get current rounding mode
inline int fegetround() {
unsigned short fpu_mode;
STREFLOP_FSTCW(fpu_mode);
return fpu_mode & 0x0C00;
}
/// Set a new rounding mode
inline int fesetround(FPU_RoundMode roundMode) {
unsigned short fpu_mode;
STREFLOP_FSTCW(fpu_mode);
fpu_mode &= 0xF3FF; // clear current mode
fpu_mode |= roundMode; // sets new mode
STREFLOP_FLDCW(fpu_mode);
return 0;
}
typedef short int fenv_t;
/// Default env. Defined in Math.cpp to be 0, and initalized on first use to the permanent holder
extern fenv_t FE_DFL_ENV;
/// Get FP env into the given structure
inline int fegetenv(fenv_t *envp) {
// check that default env exists, otherwise save it now
if (!FE_DFL_ENV) STREFLOP_FSTCW(FE_DFL_ENV);
// Now store env into argument
STREFLOP_FSTCW(*envp);
return 0;
}
/// Sets FP env from the given structure
inline int fesetenv(const fenv_t *envp) {
// check that default env exists, otherwise save it now
if (!FE_DFL_ENV) STREFLOP_FSTCW(FE_DFL_ENV);
// Now overwrite current env by argument
STREFLOP_FLDCW(*envp);
return 0;
}
/// get env and clear exceptions
inline int feholdexcept(fenv_t *envp) {
fegetenv(envp);
feclearexcept(FE_ALL_EXCEPT);
return 0;
}
template<typename T> inline void streflop_init() {
struct X {};
X Unknown_numeric_type;
// unknown types do not compile
T error = Unknown_numeric_type;
}
/// Initialize the FPU for the different types
/// this may also be called to switch between code sections using
/// different precisions
template<> inline void streflop_init<Simple>() {
unsigned short fpu_mode;
STREFLOP_FSTCW(fpu_mode);
fpu_mode &= 0xFCFF; // 32 bits internal operations
STREFLOP_FLDCW(fpu_mode);
// Enable signaling nans if compiled with this option.
#if defined(__SUPPORT_SNAN__) && !defined(USE_GML)
feraiseexcept(streflop::FPU_Exceptions(FE_INVALID | FE_DIVBYZERO | FE_OVERFLOW));
#endif
}
template<> inline void streflop_init<Double>() {
unsigned short fpu_mode;
STREFLOP_FSTCW(fpu_mode);
fpu_mode &= 0xFCFF;
fpu_mode |= 0x0200; // 64 bits internal operations
STREFLOP_FLDCW(fpu_mode);
#if defined(__SUPPORT_SNAN__) && !defined(USE_GML)
feraiseexcept(streflop::FPU_Exceptions(FE_INVALID | FE_DIVBYZERO | FE_OVERFLOW));
#endif
}
#if defined(Extended)
template<> inline void streflop_init<Extended>() {
unsigned short fpu_mode;
STREFLOP_FSTCW(fpu_mode);
fpu_mode &= 0xFCFF;
fpu_mode |= 0x0300; // 80 bits internal operations
STREFLOP_FLDCW(fpu_mode);
#if defined(__SUPPORT_SNAN__) && !defined(USE_GML)
feraiseexcept(streflop::FPU_Exceptions(FE_INVALID | FE_DIVBYZERO | FE_OVERFLOW));
#endif
}
#endif // defined(Extended)
#elif defined(STREFLOP_SSE)
/// Raise exception for these flags
inline int feraiseexcept(FPU_Exceptions excepts) {
// Just in case the compiler would store a value on the st(x) registers
unsigned short x87_mode;
STREFLOP_FSTCW(x87_mode);
x87_mode &= ~( excepts ); // generate error for selection
STREFLOP_FLDCW(x87_mode);
int sse_mode;
STREFLOP_STMXCSR(sse_mode);
sse_mode &= ~( excepts << 7 ); // generate error for selection
STREFLOP_LDMXCSR(sse_mode);
return 0;
}
/// Clear exceptions for these flags
inline int feclearexcept(int excepts) {
// Just in case the compiler would store a value on the st(x) registers
unsigned short x87_mode;
STREFLOP_FSTCW(x87_mode);
x87_mode |= excepts;
STREFLOP_FLDCW(x87_mode);
int sse_mode;
STREFLOP_STMXCSR(sse_mode);
sse_mode |= excepts << 7;
STREFLOP_LDMXCSR(sse_mode);
return 0;
}
/// Get current rounding mode
inline int fegetround() {
int sse_mode;
STREFLOP_STMXCSR(sse_mode);
return (sse_mode>>3) & 0x00000C00;
}
/// Set a new rounding mode
inline int fesetround(FPU_RoundMode roundMode) {
int sse_mode;
STREFLOP_STMXCSR(sse_mode);
sse_mode &= 0xFFFF9FFF; // clear current mode
sse_mode |= roundMode<<3; // sets new mode
STREFLOP_LDMXCSR(sse_mode);
return 0;
}
/// stores both x87 and SSE words
struct fenv_t {
int sse_mode;
short int x87_mode;
};
/// Default env. Defined in Math.cpp, structs are initialized to 0
extern fenv_t FE_DFL_ENV;
/// Get FP env into the given structure
inline int fegetenv(fenv_t *envp) {
// check that default env exists, otherwise save it now
if (!FE_DFL_ENV.x87_mode) STREFLOP_FSTCW(FE_DFL_ENV.x87_mode);
// Now store env into argument
STREFLOP_FSTCW(envp->x87_mode);
// For SSE
if (!FE_DFL_ENV.sse_mode) STREFLOP_STMXCSR(FE_DFL_ENV.sse_mode);
// Now store env into argument
STREFLOP_STMXCSR(envp->sse_mode);
return 0;
}
/// Sets FP env from the given structure
inline int fesetenv(const fenv_t *envp) {
// check that default env exists, otherwise save it now
if (!FE_DFL_ENV.x87_mode) STREFLOP_FSTCW(FE_DFL_ENV.x87_mode);
// Now overwrite current env by argument
STREFLOP_FLDCW(envp->x87_mode);
// For SSE
if (!FE_DFL_ENV.sse_mode) STREFLOP_STMXCSR(FE_DFL_ENV.sse_mode);
// Now overwrite current env by argument
STREFLOP_LDMXCSR(envp->sse_mode);
return 0;
}
/// get env and clear exceptions
inline int feholdexcept(fenv_t *envp) {
fegetenv(envp);
feclearexcept(FE_ALL_EXCEPT);
return 0;
}
template<typename T> inline void streflop_init() {
// Do nothing by default, or for unknown types
}
/// Initialize the FPU for the different types
/// this may also be called to switch between code sections using
/// different precisions
template<> inline void streflop_init<Simple>() {
// Just in case the compiler would store a value on the st(x) registers
unsigned short x87_mode;
STREFLOP_FSTCW(x87_mode);
x87_mode &= 0xFCFF; // 32 bits internal operations
STREFLOP_FLDCW(x87_mode);
int sse_mode;
STREFLOP_STMXCSR(sse_mode);
#if defined(STREFLOP_NO_DENORMALS)
sse_mode |= 0x8040; // set DAZ and FTZ
#else
sse_mode &= 0xFFFF7FBF; // clear DAZ and FTZ
#endif
STREFLOP_LDMXCSR(sse_mode);
}
template<> inline void streflop_init<Double>() {
// Just in case the compiler would store a value on the st(x) registers
unsigned short x87_mode;
STREFLOP_FSTCW(x87_mode);
x87_mode &= 0xFCFF;
x87_mode |= 0x0200; // 64 bits internal operations
STREFLOP_FLDCW(x87_mode);
int sse_mode;
STREFLOP_STMXCSR(sse_mode);
#if defined(STREFLOP_NO_DENORMALS)
sse_mode |= 0x8040; // set DAZ and FTZ
#else
sse_mode &= 0xFFFF7FBF; // clear DAZ and FTZ
#endif
STREFLOP_LDMXCSR(sse_mode);
}
#if defined(Extended)
template<> inline void streflop_init<Extended>() {
// Just in case the compiler would store a value on the st(x) registers
unsigned short x87_mode;
STREFLOP_FSTCW(x87_mode);
x87_mode &= 0xFCFF;
x87_mode |= 0x0300; // 80 bits internal operations
STREFLOP_FLDCW(x87_mode);
int sse_mode;
STREFLOP_STMXCSR(sse_mode);
#if defined(STREFLOP_NO_DENORMALS)
sse_mode |= 0x8040; // set DAZ and FTZ
#else
sse_mode &= 0xFFFF7FBF; // clear DAZ and FTZ
#endif
STREFLOP_LDMXCSR(sse_mode);
}
#endif // defined(Extended)
#elif defined(STREFLOP_SOFT)
/// Raise exception for these flags
inline int feraiseexcept(FPU_Exceptions excepts) {
// Use positive logic
SoftFloat::float_exception_realtraps |= excepts;
return 0;
}
/// Clear exceptions for these flags
inline int feclearexcept(int excepts) {
// Use positive logic
SoftFloat::float_exception_realtraps &= ~( excepts );
return 0;
}
/// Get current rounding mode
inline int fegetround() {
// see softfloat.h for the definition
switch (SoftFloat::float_rounding_mode) {
case SoftFloat::float_round_down: return FE_DOWNWARD;
case SoftFloat::float_round_up: return FE_UPWARD;
case SoftFloat::float_round_to_zero: return FE_TOWARDZERO;
default:; // is also initial mode
}
// case SoftFloat::float_round_nearest_even:
return FE_TONEAREST;
}
/// Set a new rounding mode
inline int fesetround(FPU_RoundMode roundMode) {
// see softfloat.h for the definition
switch (roundMode) {
case FE_DOWNWARD: SoftFloat::float_rounding_mode = SoftFloat::float_round_down; return 0;
case FE_UPWARD: SoftFloat::float_rounding_mode = SoftFloat::float_round_up; return 0;
case FE_TOWARDZERO: SoftFloat::float_rounding_mode = SoftFloat::float_round_to_zero; return 0;
case FE_TONEAREST: SoftFloat::float_rounding_mode = SoftFloat::float_round_nearest_even; return 0;
}
// Error, invalid mode
return 1;
}
/// SoftFloat environment comprises non-volatile state variables
struct fenv_t {
char tininess;
char rounding_mode;
int exception_realtraps;
};
/// Default env. Defined in Math.cpp, initialized to some invalid value for detection
extern fenv_t FE_DFL_ENV;
/// Get FP env into the given structure
inline int fegetenv(fenv_t *envp) {
// check that default env exists, otherwise save it now
if (FE_DFL_ENV.tininess==42) {
// First use: save default environment now
FE_DFL_ENV.tininess = SoftFloat::float_detect_tininess;
FE_DFL_ENV.rounding_mode = SoftFloat::float_rounding_mode;
FE_DFL_ENV.exception_realtraps = SoftFloat::float_exception_realtraps;
}
// Now get the current env in the given argument
envp->tininess = SoftFloat::float_detect_tininess;
envp->rounding_mode = SoftFloat::float_rounding_mode;
envp->exception_realtraps = SoftFloat::float_exception_realtraps;
return 0;
}
/// Sets FP env from the given structure
inline int fesetenv(const fenv_t *envp) {
// check that default env exists, otherwise save it now
if (FE_DFL_ENV.tininess==42) {
// First use: save default environment now
FE_DFL_ENV.tininess = SoftFloat::float_detect_tininess;
FE_DFL_ENV.rounding_mode = SoftFloat::float_rounding_mode;
FE_DFL_ENV.exception_realtraps = SoftFloat::float_exception_realtraps;
}
// Now get the current env in the given argument
SoftFloat::float_detect_tininess = envp->tininess;
SoftFloat::float_rounding_mode = envp->rounding_mode;
SoftFloat::float_exception_realtraps = envp->exception_realtraps;
return 0;
}
/// get env and clear exceptions
inline int feholdexcept(fenv_t *envp) {
fegetenv(envp);
feclearexcept(FE_ALL_EXCEPT);
return 0;
}
template<typename T> inline void streflop_init() {
// Do nothing by default, or for unknown types
}
/// Initialize the FPU for the different types
/// this may also be called to switch between code sections using
/// different precisions
template<> inline void streflop_init<Simple>() {
}
template<> inline void streflop_init<Double>() {
}
template<> inline void streflop_init<Extended>() {
}
#else // defined(STREFLOP_X87)
#error STREFLOP: Invalid combination or unknown FPU type.
#endif // defined(STREFLOP_X87)
}
#endif // STREFLOP_FPU_H

139
source/shared_lib/include/streflop/IntegerTypes.h Normal file
View File

@@ -0,0 +1,139 @@
/*
streflop: STandalone REproducible FLOating-Point
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
Heavily relies on GNU Libm, itself depending on netlib fplibm, GNU MP, and IBM MP lib.
Uses SoftFloat too.
Please read the history and copyright information in the documentation provided with the source code
*/
#ifndef STREFLOP_INTEGER_TYPES_H
#define STREFLOP_INTEGER_TYPES_H
// Given that:
// - C++ template metaprogramming is Turing complete
// - types are variables and integers are constants
// - sizeof(type) is known by the C++ compiler and a constant
// - sizeof(char) is 1 by definition, even if char != 8 bits
// Then: It is possible to derive the sized ints at compile time in C++, unlike C
// Note: This is NOT the same as the int32_t, etc from C99, in case char != 8 bits
// For this reason, redefine the macro below if this is not the case for you
// Note2: Even if char != 8 bits, it's still possible to define ints in terms of number of char!
#define STREFLOP_INTEGER_TYPES_CHAR_BITS 8
// Avoid conflict with system types, if any
namespace streflop {
// Template meta-programming: this is the "program" which variables are types and constants are int template arguments
// Algorithm: provide an expected size, and recursively increase the integer types till the size match
template<typename int_type, int expected_size, bool final_recursion> struct SizedTypeMaker {
};
// start by long long to provide the recursion terminal condition
// false : the expected_size does not exist: do not define a type, there will be a compilation error
template<int expected_size> struct SizedTypeMaker<long long, expected_size, false> {
// Error: Integer type with expected size does not exist
};
// true: end recursion by defining the correct type to be long long
template<int expected_size> struct SizedTypeMaker<long long, expected_size, true> {
typedef long long final_type;
};
// false : recurse by increasing the integer type till it reaches the expected size
template<int expected_size> struct SizedTypeMaker<long, expected_size, false> {
typedef typename SizedTypeMaker<long long, expected_size, (sizeof(long long)==expected_size)>::final_type final_type;
};
// true: end recursion by defining the correct type to be long
template<int expected_size> struct SizedTypeMaker<long, expected_size, true> {
typedef long final_type;
};
// false : recurse by increasing the integer type till it reaches the expected size
template<int expected_size> struct SizedTypeMaker<int, expected_size, false> {
typedef typename SizedTypeMaker<long, expected_size, (sizeof(long)==expected_size)>::final_type final_type;
};
// true: end recursion by defining the correct type to be int
template<int expected_size> struct SizedTypeMaker<int, expected_size, true> {
typedef int final_type;
};
// false : recurse by increasing the integer type till it reaches the expected size
template<int expected_size> struct SizedTypeMaker<short, expected_size, false> {
typedef typename SizedTypeMaker<int, expected_size, (sizeof(int)==expected_size)>::final_type final_type;
};
// true: end recursion by defining the correct type to be short
template<int expected_size> struct SizedTypeMaker<short, expected_size, true> {
typedef short final_type;
};
// Do it again for unsigned types
// false : the expected_size does not exist: do not define a type, there will be a compilation error
template<int expected_size> struct SizedTypeMaker<unsigned long long, expected_size, false> {
// Error: Integer type with expected size does not exist
};
// true: end recursion by defining the correct type to be long long
template<int expected_size> struct SizedTypeMaker<unsigned long long, expected_size, true> {
typedef unsigned long long final_type;
};
// false : recurse by increasing the integer type till it reaches the expected size
template<int expected_size> struct SizedTypeMaker<unsigned long, expected_size, false> {
typedef typename SizedTypeMaker<unsigned long long, expected_size, (sizeof(unsigned long long)==expected_size)>::final_type final_type;
};
// true: end recursion by defining the correct type to be long
template<int expected_size> struct SizedTypeMaker<unsigned long, expected_size, true> {
typedef unsigned long final_type;
};
// false : recurse by increasing the integer type till it reaches the expected size
template<int expected_size> struct SizedTypeMaker<unsigned int, expected_size, false> {
typedef typename SizedTypeMaker<unsigned long, expected_size, (sizeof(unsigned long)==expected_size)>::final_type final_type;
};
// true: end recursion by defining the correct type to be int
template<int expected_size> struct SizedTypeMaker<unsigned int, expected_size, true> {
typedef unsigned int final_type;
};
// false : recurse by increasing the integer type till it reaches the expected size
template<int expected_size> struct SizedTypeMaker<unsigned short, expected_size, false> {
typedef typename SizedTypeMaker<unsigned int, expected_size, (sizeof(unsigned int)==expected_size)>::final_type final_type;
};
// true: end recursion by defining the correct type to be short
template<int expected_size> struct SizedTypeMaker<unsigned short, expected_size, true> {
typedef unsigned short final_type;
};
// Utility to get an int type with the selected size IN BITS
template<int N> struct SizedInteger {
typedef typename SizedTypeMaker<short, (N/STREFLOP_INTEGER_TYPES_CHAR_BITS), (sizeof(short)==(N/STREFLOP_INTEGER_TYPES_CHAR_BITS))>::final_type Type;
};
template<int N> struct SizedUnsignedInteger {
typedef typename SizedTypeMaker<unsigned short, (N/STREFLOP_INTEGER_TYPES_CHAR_BITS), (sizeof(unsigned short)==(N/STREFLOP_INTEGER_TYPES_CHAR_BITS))>::final_type Type;
};
// Specialize for size = STREFLOP_INTEGER_TYPES_CHAR_BITS
template<> struct SizedInteger<STREFLOP_INTEGER_TYPES_CHAR_BITS> {
typedef char Type;
};
template<> struct SizedUnsignedInteger<STREFLOP_INTEGER_TYPES_CHAR_BITS> {
typedef unsigned char Type;
};
}
#endif

510
source/shared_lib/include/streflop/LGPL.txt Normal file
View File

@@ -0,0 +1,510 @@
GNU LESSER GENERAL PUBLIC LICENSE
Version 2.1, February 1999
Copyright (C) 1991, 1999 Free Software Foundation, Inc.
59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
Everyone is permitted to copy and distribute verbatim copies
of this license document, but changing it is not allowed.
[This is the first released version of the Lesser GPL. It also counts
as the successor of the GNU Library Public License, version 2, hence
the version number 2.1.]
Preamble
The licenses for most software are designed to take away your
freedom to share and change it. By contrast, the GNU General Public
Licenses are intended to guarantee your freedom to share and change
free software--to make sure the software is free for all its users.
This license, the Lesser General Public License, applies to some
specially designated software packages--typically libraries--of the
Free Software Foundation and other authors who decide to use it. You
can use it too, but we suggest you first think carefully about whether
this license or the ordinary General Public License is the better
strategy to use in any particular case, based on the explanations
below.
When we speak of free software, we are referring to freedom of use,
not price. Our General Public Licenses are designed to make sure that
you have the freedom to distribute copies of free software (and charge
for this service if you wish); that you receive source code or can get
it if you want it; that you can change the software and use pieces of
it in new free programs; and that you are informed that you can do
these things.
To protect your rights, we need to make restrictions that forbid
distributors to deny you these rights or to ask you to surrender these
rights. These restrictions translate to certain responsibilities for
you if you distribute copies of the library or if you modify it.
For example, if you distribute copies of the library, whether gratis
or for a fee, you must give the recipients all the rights that we gave
you. You must make sure that they, too, receive or can get the source
code. If you link other code with the library, you must provide
complete object files to the recipients, so that they can relink them
with the library after making changes to the library and recompiling
it. And you must show them these terms so they know their rights.
We protect your rights with a two-step method: (1) we copyright the
library, and (2) we offer you this license, which gives you legal
permission to copy, distribute and/or modify the library.
To protect each distributor, we want to make it very clear that
there is no warranty for the free library. Also, if the library is
modified by someone else and passed on, the recipients should know
that what they have is not the original version, so that the original
author's reputation will not be affected by problems that might be
introduced by others.
^L
Finally, software patents pose a constant threat to the existence of
any free program. We wish to make sure that a company cannot
effectively restrict the users of a free program by obtaining a
restrictive license from a patent holder. Therefore, we insist that
any patent license obtained for a version of the library must be
consistent with the full freedom of use specified in this license.
Most GNU software, including some libraries, is covered by the
ordinary GNU General Public License. This license, the GNU Lesser
General Public License, applies to certain designated libraries, and
is quite different from the ordinary General Public License. We use
this license for certain libraries in order to permit linking those
libraries into non-free programs.
When a program is linked with a library, whether statically or using
a shared library, the combination of the two is legally speaking a
combined work, a derivative of the original library. The ordinary
General Public License therefore permits such linking only if the
entire combination fits its criteria of freedom. The Lesser General
Public License permits more lax criteria for linking other code with
the library.
We call this license the "Lesser" General Public License because it
does Less to protect the user's freedom than the ordinary General
Public License. It also provides other free software developers Less
of an advantage over competing non-free programs. These disadvantages
are the reason we use the ordinary General Public License for many
libraries. However, the Lesser license provides advantages in certain
special circumstances.
For example, on rare occasions, there may be a special need to
encourage the widest possible use of a certain library, so that it
becomes a de-facto standard. To achieve this, non-free programs must
be allowed to use the library. A more frequent case is that a free
library does the same job as widely used non-free libraries. In this
case, there is little to gain by limiting the free library to free
software only, so we use the Lesser General Public License.
In other cases, permission to use a particular library in non-free
programs enables a greater number of people to use a large body of
free software. For example, permission to use the GNU C Library in
non-free programs enables many more people to use the whole GNU
operating system, as well as its variant, the GNU/Linux operating
system.
Although the Lesser General Public License is Less protective of the
users' freedom, it does ensure that the user of a program that is
linked with the Library has the freedom and the wherewithal to run
that program using a modified version of the Library.
The precise terms and conditions for copying, distribution and
modification follow. Pay close attention to the difference between a
"work based on the library" and a "work that uses the library". The
former contains code derived from the library, whereas the latter must
be combined with the library in order to run.
^L
GNU LESSER GENERAL PUBLIC LICENSE
TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION
0. This License Agreement applies to any software library or other
program which contains a notice placed by the copyright holder or
other authorized party saying it may be distributed under the terms of
this Lesser General Public License (also called "this License").
Each licensee is addressed as "you".
A "library" means a collection of software functions and/or data
prepared so as to be conveniently linked with application programs
(which use some of those functions and data) to form executables.
The "Library", below, refers to any such software library or work
which has been distributed under these terms. A "work based on the
Library" means either the Library or any derivative work under
copyright law: that is to say, a work containing the Library or a
portion of it, either verbatim or with modifications and/or translated
straightforwardly into another language. (Hereinafter, translation is
included without limitation in the term "modification".)
"Source code" for a work means the preferred form of the work for
making modifications to it. For a library, complete source code means
all the source code for all modules it contains, plus any associated
interface definition files, plus the scripts used to control
compilation and installation of the library.
Activities other than copying, distribution and modification are not
covered by this License; they are outside its scope. The act of
running a program using the Library is not restricted, and output from
such a program is covered only if its contents constitute a work based
on the Library (independent of the use of the Library in a tool for
writing it). Whether that is true depends on what the Library does
and what the program that uses the Library does.
1. You may copy and distribute verbatim copies of the Library's
complete source code as you receive it, in any medium, provided that
you conspicuously and appropriately publish on each copy an
appropriate copyright notice and disclaimer of warranty; keep intact
all the notices that refer to this License and to the absence of any
warranty; and distribute a copy of this License along with the
Library.
You may charge a fee for the physical act of transferring a copy,
and you may at your option offer warranty protection in exchange for a
fee.
2. You may modify your copy or copies of the Library or any portion
of it, thus forming a work based on the Library, and copy and
distribute such modifications or work under the terms of Section 1
above, provided that you also meet all of these conditions:
a) The modified work must itself be a software library.
b) You must cause the files modified to carry prominent notices
stating that you changed the files and the date of any change.
c) You must cause the whole of the work to be licensed at no
charge to all third parties under the terms of this License.
d) If a facility in the modified Library refers to a function or a
table of data to be supplied by an application program that uses
the facility, other than as an argument passed when the facility
is invoked, then you must make a good faith effort to ensure that,
in the event an application does not supply such function or
table, the facility still operates, and performs whatever part of
its purpose remains meaningful.
(For example, a function in a library to compute square roots has
a purpose that is entirely well-defined independent of the
application. Therefore, Subsection 2d requires that any
application-supplied function or table used by this function must
be optional: if the application does not supply it, the square
root function must still compute square roots.)
These requirements apply to the modified work as a whole. If
identifiable sections of that work are not derived from the Library,
and can be reasonably considered independent and separate works in
themselves, then this License, and its terms, do not apply to those
sections when you distribute them as separate works. But when you
distribute the same sections as part of a whole which is a work based
on the Library, the distribution of the whole must be on the terms of
this License, whose permissions for other licensees extend to the
entire whole, and thus to each and every part regardless of who wrote
it.
Thus, it is not the intent of this section to claim rights or contest
your rights to work written entirely by you; rather, the intent is to
exercise the right to control the distribution of derivative or
collective works based on the Library.
In addition, mere aggregation of another work not based on the Library
with the Library (or with a work based on the Library) on a volume of
a storage or distribution medium does not bring the other work under
the scope of this License.
3. You may opt to apply the terms of the ordinary GNU General Public
License instead of this License to a given copy of the Library. To do
this, you must alter all the notices that refer to this License, so
that they refer to the ordinary GNU General Public License, version 2,
instead of to this License. (If a newer version than version 2 of the
ordinary GNU General Public License has appeared, then you can specify
that version instead if you wish.) Do not make any other change in
these notices.
^L
Once this change is made in a given copy, it is irreversible for
that copy, so the ordinary GNU General Public License applies to all
subsequent copies and derivative works made from that copy.
This option is useful when you wish to copy part of the code of
the Library into a program that is not a library.
4. You may copy and distribute the Library (or a portion or
derivative of it, under Section 2) in object code or executable form
under the terms of Sections 1 and 2 above provided that you accompany
it with the complete corresponding machine-readable source code, which
must be distributed under the terms of Sections 1 and 2 above on a
medium customarily used for software interchange.
If distribution of object code is made by offering access to copy
from a designated place, then offering equivalent access to copy the
source code from the same place satisfies the requirement to
distribute the source code, even though third parties are not
compelled to copy the source along with the object code.
5. A program that contains no derivative of any portion of the
Library, but is designed to work with the Library by being compiled or
linked with it, is called a "work that uses the Library". Such a
work, in isolation, is not a derivative work of the Library, and
therefore falls outside the scope of this License.
However, linking a "work that uses the Library" with the Library
creates an executable that is a derivative of the Library (because it
contains portions of the Library), rather than a "work that uses the
library". The executable is therefore covered by this License.
Section 6 states terms for distribution of such executables.
When a "work that uses the Library" uses material from a header file
that is part of the Library, the object code for the work may be a
derivative work of the Library even though the source code is not.
Whether this is true is especially significant if the work can be
linked without the Library, or if the work is itself a library. The
threshold for this to be true is not precisely defined by law.
If such an object file uses only numerical parameters, data
structure layouts and accessors, and small macros and small inline
functions (ten lines or less in length), then the use of the object
file is unrestricted, regardless of whether it is legally a derivative
work. (Executables containing this object code plus portions of the
Library will still fall under Section 6.)
Otherwise, if the work is a derivative of the Library, you may
distribute the object code for the work under the terms of Section 6.
Any executables containing that work also fall under Section 6,
whether or not they are linked directly with the Library itself.
^L
6. As an exception to the Sections above, you may also combine or
link a "work that uses the Library" with the Library to produce a
work containing portions of the Library, and distribute that work
under terms of your choice, provided that the terms permit
modification of the work for the customer's own use and reverse
engineering for debugging such modifications.
You must give prominent notice with each copy of the work that the
Library is used in it and that the Library and its use are covered by
this License. You must supply a copy of this License. If the work
during execution displays copyright notices, you must include the
copyright notice for the Library among them, as well as a reference
directing the user to the copy of this License. Also, you must do one
of these things:
a) Accompany the work with the complete corresponding
machine-readable source code for the Library including whatever
changes were used in the work (which must be distributed under
Sections 1 and 2 above); and, if the work is an executable linked
with the Library, with the complete machine-readable "work that
uses the Library", as object code and/or source code, so that the
user can modify the Library and then relink to produce a modified
executable containing the modified Library. (It is understood
that the user who changes the contents of definitions files in the
Library will not necessarily be able to recompile the application
to use the modified definitions.)
b) Use a suitable shared library mechanism for linking with the
Library. A suitable mechanism is one that (1) uses at run time a
copy of the library already present on the user's computer system,
rather than copying library functions into the executable, and (2)
will operate properly with a modified version of the library, if
the user installs one, as long as the modified version is
interface-compatible with the version that the work was made with.
c) Accompany the work with a written offer, valid for at least
three years, to give the same user the materials specified in
Subsection 6a, above, for a charge no more than the cost of
performing this distribution.
d) If distribution of the work is made by offering access to copy
from a designated place, offer equivalent access to copy the above
specified materials from the same place.
e) Verify that the user has already received a copy of these
materials or that you have already sent this user a copy.
For an executable, the required form of the "work that uses the
Library" must include any data and utility programs needed for
reproducing the executable from it. However, as a special exception,
the materials to be distributed need not include anything that is
normally distributed (in either source or binary form) with the major
components (compiler, kernel, and so on) of the operating system on
which the executable runs, unless that component itself accompanies
the executable.
It may happen that this requirement contradicts the license
restrictions of other proprietary libraries that do not normally
accompany the operating system. Such a contradiction means you cannot
use both them and the Library together in an executable that you
distribute.
^L
7. You may place library facilities that are a work based on the
Library side-by-side in a single library together with other library
facilities not covered by this License, and distribute such a combined
library, provided that the separate distribution of the work based on
the Library and of the other library facilities is otherwise
permitted, and provided that you do these two things:
a) Accompany the combined library with a copy of the same work
based on the Library, uncombined with any other library
facilities. This must be distributed under the terms of the
Sections above.
b) Give prominent notice with the combined library of the fact
that part of it is a work based on the Library, and explaining
where to find the accompanying uncombined form of the same work.
8. You may not copy, modify, sublicense, link with, or distribute
the Library except as expressly provided under this License. Any
attempt otherwise to copy, modify, sublicense, link with, or
distribute the Library is void, and will automatically terminate your
rights under this License. However, parties who have received copies,
or rights, from you under this License will not have their licenses
terminated so long as such parties remain in full compliance.
9. You are not required to accept this License, since you have not
signed it. However, nothing else grants you permission to modify or
distribute the Library or its derivative works. These actions are
prohibited by law if you do not accept this License. Therefore, by
modifying or distributing the Library (or any work based on the
Library), you indicate your acceptance of this License to do so, and
all its terms and conditions for copying, distributing or modifying
the Library or works based on it.
10. Each time you redistribute the Library (or any work based on the
Library), the recipient automatically receives a license from the
original licensor to copy, distribute, link with or modify the Library
subject to these terms and conditions. You may not impose any further
restrictions on the recipients' exercise of the rights granted herein.
You are not responsible for enforcing compliance by third parties with
this License.
^L
11. If, as a consequence of a court judgment or allegation of patent
infringement or for any other reason (not limited to patent issues),
conditions are imposed on you (whether by court order, agreement or
otherwise) that contradict the conditions of this License, they do not
excuse you from the conditions of this License. If you cannot
distribute so as to satisfy simultaneously your obligations under this
License and any other pertinent obligations, then as a consequence you
may not distribute the Library at all. For example, if a patent
license would not permit royalty-free redistribution of the Library by
all those who receive copies directly or indirectly through you, then
the only way you could satisfy both it and this License would be to
refrain entirely from distribution of the Library.
If any portion of this section is held invalid or unenforceable under
any particular circumstance, the balance of the section is intended to
apply, and the section as a whole is intended to apply in other
circumstances.
It is not the purpose of this section to induce you to infringe any
patents or other property right claims or to contest validity of any
such claims; this section has the sole purpose of protecting the
integrity of the free software distribution system which is
implemented by public license practices. Many people have made
generous contributions to the wide range of software distributed
through that system in reliance on consistent application of that
system; it is up to the author/donor to decide if he or she is willing
to distribute software through any other system and a licensee cannot
impose that choice.
This section is intended to make thoroughly clear what is believed to
be a consequence of the rest of this License.
12. If the distribution and/or use of the Library is restricted in
certain countries either by patents or by copyrighted interfaces, the
original copyright holder who places the Library under this License
may add an explicit geographical distribution limitation excluding those
countries, so that distribution is permitted only in or among
countries not thus excluded. In such case, this License incorporates
the limitation as if written in the body of this License.
13. The Free Software Foundation may publish revised and/or new
versions of the Lesser General Public License from time to time.
Such new versions will be similar in spirit to the present version,
but may differ in detail to address new problems or concerns.
Each version is given a distinguishing version number. If the Library
specifies a version number of this License which applies to it and
"any later version", you have the option of following the terms and
conditions either of that version or of any later version published by
the Free Software Foundation. If the Library does not specify a
license version number, you may choose any version ever published by
the Free Software Foundation.
^L
14. If you wish to incorporate parts of the Library into other free
programs whose distribution conditions are incompatible with these,
write to the author to ask for permission. For software which is
copyrighted by the Free Software Foundation, write to the Free
Software Foundation; we sometimes make exceptions for this. Our
decision will be guided by the two goals of preserving the free status
of all derivatives of our free software and of promoting the sharing
and reuse of software generally.
NO WARRANTY
15. BECAUSE THE LIBRARY IS LICENSED FREE OF CHARGE, THERE IS NO
WARRANTY FOR THE LIBRARY, TO THE EXTENT PERMITTED BY APPLICABLE LAW.
EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR
OTHER PARTIES PROVIDE THE LIBRARY "AS IS" WITHOUT WARRANTY OF ANY
KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE
LIBRARY IS WITH YOU. SHOULD THE LIBRARY PROVE DEFECTIVE, YOU ASSUME
THE COST OF ALL NECESSARY SERVICING, REPAIR OR CORRECTION.
16. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN
WRITING WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY
AND/OR REDISTRIBUTE THE LIBRARY AS PERMITTED ABOVE, BE LIABLE TO YOU
FOR DAMAGES, INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR
CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OR INABILITY TO USE THE
LIBRARY (INCLUDING BUT NOT LIMITED TO LOSS OF DATA OR DATA BEING
RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD PARTIES OR A
FAILURE OF THE LIBRARY TO OPERATE WITH ANY OTHER SOFTWARE), EVEN IF
SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH
DAMAGES.
END OF TERMS AND CONDITIONS
^L
How to Apply These Terms to Your New Libraries
If you develop a new library, and you want it to be of the greatest
possible use to the public, we recommend making it free software that
everyone can redistribute and change. You can do so by permitting
redistribution under these terms (or, alternatively, under the terms
of the ordinary General Public License).
To apply these terms, attach the following notices to the library.
It is safest to attach them to the start of each source file to most
effectively convey the exclusion of warranty; and each file should
have at least the "copyright" line and a pointer to where the full
notice is found.
<one line to give the library's name and a brief idea of what it does.>
Copyright (C) <year> <name of author>
This library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with this library; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
Also add information on how to contact you by electronic and paper mail.
You should also get your employer (if you work as a programmer) or
your school, if any, to sign a "copyright disclaimer" for the library,
if necessary. Here is a sample; alter the names:
Yoyodyne, Inc., hereby disclaims all copyright interest in the
library `Frob' (a library for tweaking knobs) written by James
Random Hacker.
<signature of Ty Coon>, 1 April 1990
Ty Coon, President of Vice
That's all there is to it!

205
source/shared_lib/include/streflop/README.txt Normal file
View File

@@ -0,0 +1,205 @@
STandalone REproducible FLOating-Point library
Version 0.2, june 2006.
Nicolas Brodu. See also the acknowledgments below.
For a quick setup guide, see the "usage" sections below.
Presentation:
Floating-point computations are strongly dependent on the FPU hardware implementation, the compiler and its optimizations, and the system mathematical library (libm). Experiments are usually reproducible only on the same machine with the same system library and the same compiler using the same options. Even then, the C++ standard does NOT guarantee reproducible results. Example:
double x = 1.0; x /= 10.0;
double y = x;
...
if (y == x) {
// THIS IS NOT ALWAYS TRUE!
}
A related but more general problem is random number generation, often necessary for experiment reproducibility:
set_random_seed(42);
double x = random_number();
May not return the same x across different FPUs / systems / compilers / etc.
These problems are related to:
- Some FPU (like x87 on Linux) keep by default an internal precision (80 bits) larger than the type size in memory. In the first example, if x is on the stack but y in a register, the comparison will fail. Worse, whether x and/or y is on the stack or register depends on what you put in the ... section. This problem can be solved by restricting the internal FPU precision to the type size in memory. Unfortunately, on x87, some operations like the transcendental functions are always computed on 80 bits...
- How well your FPU implements the IEEE754 standard, and in particular, denormal number operations. The provided arithmetic_test code reproducibly gives different results between SSE and x87. This problem is NOT generally solved by restricting the internal precision to the type size in memory. There may be relations however, in particular a denormal float may be a normal double, so the precision matters.
- How your compiler interacts with your system math library (libm) especially with optimization flags. In particular, gcc-3.4 and gcc-4.0 series may give different results on the same system with the same library. This problem is partially solved by changing the compiler options.
- Even then, the IEEE754 standard has some loopholes concerning NaN types. In particular, when serializing results, binary files may differ. This problem is solved by comparing numerical values and not bit patterns.
More points to consider:
- SSE has an option to accelerate computations by approximating denormals by 0. For some applications this is plain wrong, but for other applications like real-time DSP, denormals are a real pain and this is just what's needed (reports have been made that a denormal multiply may take as much as 30 times a normal multiply). Due to the default round-to-nearest mode, the denormals tend NOT to cancel to 0 and may instead accumulate. Deactivating denormals is built-in the core SSE FPU, but unfortunately that's not reproducible on x87. In that case, it's possible to check for denormal conditions after each operation (including assignment) and then flush to zero, thanks to a wrapper type, for reproducibility (but at the expense of performance on x87).
- No external dependency. The more dependencies, the more chances of a version mismatch. This library should be standalone, providing the whole libm features without resorting to system specific includes or other packages. This way, it can be included in a project as is, with minimal specialization (and risk of misconfiguration).
Proposed solution:
- Provide Simple, Double and Extended types that match the native types (float, double, long double), but that additionally take care of the FPU internal precision and denormal handling. These types may be simple aliases or C++ wrappers, depending on the FPU and configuration. Note: Extended is the 80-bit type defined by x87.
- Reimplement the libm as a standalone code that uses these types. Note: For extended support, see below.
- As a bonus, provide a random number generator so as to make standalone experiments reproducible.
- Compare with a software floating-point reference implementation.
Usage (programming):
- Include "streflop.h" in place of <math.h>, and link with streflop.a instead of libm. All streflop functions are protected by a namespace, so if another part of your program uses libm there is no risk of confusion at link time. However, including the correct file matters.
- Use the streflop namespace, and the Simple, Double and Extended types as needed, instead of float, double, long double. The streflop types may actually be aliases to the C++ types, or wrapper classes that redefine the operators transparently.
- You should also call streflop_init<FloatType> with FloatType=Simple,Double,Extended before using that type. You should use only this type (ex: Simple) until the next call to streflop_init. That is, separate your code in blocks using one type at a time. In the simplest case, use streflop_init for your chosen type at the beginning of your program and stick to that type later on. These init functions are necessary to set the correct FPU flags. See also the notes below.
- You may have a look at arithmeticTest.cpp and randomTest.cpp for examples.
Usage (standalone build):
- Edit Makefile.common to configure which FPU/denormal setup you choose, by defining one of STREFLOP_SSE, STREFLOP_X87, STREFLOP_SOFT, and optionally STREFLOP_NO_DENORMALS. See the configurations grid below.
- If you're using the software floating-point implementation on a big-endian machine, change the System.h file accordingly. If your target system size has a char type larger than 8 bits, then check Integer.h. In both cases you're on your own (this is untested).
- Check the notes below before changing the compiler options.
Usage (including in a project):
- Copy the whole streflop source somewhere in your project, for example as a "streflop" subdirectory.
- Check the steps for the standalone build. Potentially automate the choices using your build system, like autoconf.
- Call "make -C streflop" somewhere in you own build process, or integrate the streflop build in your build system.
- See the programming usage above. Don't forget to "-Istreflop" and "-Lstreflop -lstreflop.a" when compiling and linking.
Configurations grid:
| SSE | x87 | Soft |
------------+---------+-----------+----------+
denormals | Simple *| Simple *| Simple |
| Double *| Double *| Double |
| | Extended *| Extended |
------------+---------+-----------+----------+
no denormal | Simple *| Simple |
| Double *| Double |
| | Extended |
------------+---------+-----------+
One cell in this grid must be selected at configure time. All types within that cell are then available at compile and run time.
Types marked * are aliases to the native float/double/long double, with support by FPU flags.
The other types are wrapper classes that behave like the native types.
Apart for the bit representation of NaN values:
- "Denormals SSE / Denormals Soft" with the same precision should give the same results.
- "No denormal SSE / No denormal x87" with the same precision should give the same results.
- "Denormals x87 extended / Denormals Soft extended" should give the same results.
- "Denormals SSE / Denormals x87" with the same precision may differ but only for some unfrequent occurences involving denormal numbers.
- All other configurations give different results.
Comparison criteria:
- Best performance is achieved by "no denormal/SSE simple". What matters most for performance is wrapper/native, the size, then denormals or not (unless using lots of denormals, in which case "no denormals" may matter more than size or wrapper).
- Best precision is achieved by "denormals extended" modes. What matters most for precision is the the size, then denormal or not.
- Best IEEE754 conformance is achieved by "Soft" modes (and equivalent results above with SSE/x87). Conformance is achieved only for denormals.
Notes:
- Beware of too aggressive optimization options! In particular, since this code relies on reinterpret_cast and unions, the compiler must not assume strict aliasing. For g++ optimization levels 2 and 3, this assumption is unfortunately the default. Similarly, the compiler should not assume that NaN can be ignored, or that the FPU has a constant rounding mode. Ex: -O3 -fno-strict-aliasing -frounding-math -fsignaling-nans.
- You should also set correct FPU options, like -mfpmath=sse -msse -msse2. The -msse2 is important, there are cases where g++ refuses to use SSE (and silently falls back to x87) when using -msse and not -msse2. This also means you cannot reliably use this library with gcc on systems where only sse (but not sse2) is present, like some athlon-xp cores.
- The system libm will almost surely produce different numerical results depending on your FPU, compiler and options, etc. The rationale is, using this library will increase the reproducibility of your experiments compared to using the system libm. If you want guaranteed (but slower) reproducible results across all machine configurations, without caring for denormals or whatever else, then use a multiprecision software library like GNU MP. If you want to use the hardware FPU in a controlled environment that can retain some reproducibility, then use this library. If it does not fit your needs, then improve it: After all, this is free software :)
- The following C99 trap and rounding mode functions are implemented, even for the software floating-point implementation: fe(get|set)round, fe(get|set)env, and feholdexcept. You may call them to change rounding modes and to trap special conditions. These functions are expected to work correctly, insofar as the FPU works as intended*, but they have not been extensively tested.
* in particular, reports have been made that the x87 FPU denormal trap sometimes fails.
- Really beware of aggressive optimization! Separate your code into INDEPENDENT BLOCKS. I mean it. This code is wrong:
streflop_init<Simple>();
Simple s = (1.0/4294967295.0);
displayHex(cout, s) << endl;
streflop_init<Double>();
Double d = (1.0/4294967295.0);
displayHex(cout, d) << endl;
THIS NOT WORK CORRECTLY in -O3, but will do fine in -O0. This is because the compiler "optimizes" the constant computation only once for both lines in -O3, which is plain wrong since the precision is different. The only way to ensure this does not happen is to separate your code in logical units:
void func1() {
streflop_init<Simple>();
Simple s = (1.0/4294967295.0);
displayHex(cout, s) << endl;
}
void func2() {
streflop_init<Double>();
Double d = (1.0/4294967295.0);
displayHex(cout, d) << endl;
}
Even then, you'd better be careful with interprocedural optimization options. If possible, put func1 and func2 in 2 separate compilation units.
BUGS and discrepancies:
- Do what you want with this library, but at your own risks, and according to the LGPL (see the LGPL.txt file).
- There is the possibility of unknown bugs. And this is based on GNU libm 2.4, so any potential bug in that version are almost surely present in streflop too.
- Extended support is INCOMPLETE. Proper functions are missing, in particular the trigonometric functions. The ldbl-96 implementation of the libm does not contain a generic implementation for these files. Since strelop enforces strict separation of Extended and Double functions, these functions were instead implemented by temporarily switching to Double using streflop_init<Double>, calling the function and storing the result on the stack, switching back to streflop_init<Extended>, then converting the result to an Extended number.
Acknowledgments:
- This code heavily relies on GNU Libm, itself depending on Sun's netlib fplibm, GNU MP, and IBM's multi-precision library.
- This code uses the SoftFloat library for the software floating-point implementation.
- The random number generator is the Mersene Twister, created by Takuji Nishimura and Makoto Matsumoto, and adapted to C++ for this project. Please read the (BSD-like) license in Random.cpp if you intend to make binary packages of this library (and according to LGPL).
- Please read the history and copyright information in the accompanying README.txt files in the libm and softfloat directories, as well as the LGPL.txt file in this directory.
- Thanks to Tobi Vollebregt for feedback, Win32 reports, and patches.
How you can help:
- Test the library and report potential bugs.
- Port the library to new FPU and operating systems.
- Help extend the GNU libm first, and only then import that work in this project.
Nicolas Brodu, june 2006.

245
source/shared_lib/include/streflop/Random.h Normal file
View File

@@ -0,0 +1,245 @@
/*
streflop: STandalone REproducible FLOating-Point
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
Heavily relies on GNU Libm, itself depending on netlib fplibm, GNU MP, and IBM MP lib.
Uses SoftFloat too.
Please read the history and copyright information in the documentation provided with the source code
*/
#ifndef RANDOM_H
#define RANDOM_H
// Need sized integer types, which are now system-independent thanks to template metaprogramming
#include "IntegerTypes.h"
namespace streflop {
/** Random state holder object
Declare one of this per thread, and use it as context to the random functions
Object is properly set by the RandomInit functions
This allows it to remain POD type
*/
struct RandomState {
#if !defined(STREFLOP_RANDOM_GEN_SIZE)
#define STREFLOP_RANDOM_GEN_SIZE 32
#endif
// state vector
SizedUnsignedInteger<STREFLOP_RANDOM_GEN_SIZE>::Type mt[19968/STREFLOP_RANDOM_GEN_SIZE];
int mti;
// random seed that was used for initialization
SizedUnsignedInteger<32>::Type seed;
}
#ifdef __GNUC__
__attribute__ ((aligned (64))) // align state vector on cache line size
#endif
;
/// Default random state holder
extern RandomState DefaultRandomState;
/** Initialize the random number generator with the given seed.
By default, the seed is taken from system time and printed out
so the experiment is reproducible by inputing the same seed again.
You can set here a previous seed to reproduce it.
This interface allows independance from the actual RNG used,
and/or system functions.
The RNG used is the Mersenne twister implementation by the
original authors Takuji Nishimura and Makoto Matsumoto.
See also Random.cpp for more information.
*/
SizedUnsignedInteger<32>::Type RandomInit(RandomState& state = DefaultRandomState);
SizedUnsignedInteger<32>::Type RandomInit(SizedUnsignedInteger<32>::Type seed, RandomState& state = DefaultRandomState);
/// Returns the random seed that was used for the initialization
/// Defaults to 0 if the RNG is not yet initialized
SizedUnsignedInteger<32>::Type RandomSeed(RandomState& state = DefaultRandomState);
/** Returns a random number from a uniform distribution.
All integer types are supported, as well as Simple, Double, and Extended
The Random(min, max) template takes as argument:
- bool: whether or not including the min bound
- bool: whether or not including the max bound
- type: to generate a number from that type, and decide on min/max bounds
Example: Double x = Random<true, false, Double>(7.0, 18.0)
This will return a Double number between 7.0 (included) and 18.0 (excluded)
This works for both float and integer types.
Aliases are named like RandomXY with X,Y = E,I for bounds Excluded,Included.
Example: RandomEI(min,max) will return a number between min (excluded) and max (included).
The Random() template returns a number of the given type chosen uniformly between
all representable numbers of that type.
- For integer types, this means what you expect: any int, char, whatever on the whole range
- For float types, just recall that there are as many representable floats in each range
2^x - 2^(x+1), this is what floating-point means
Notes:
- If min > max then the result is undefined.
- If you ask for an empty interval (like RandomEE with min==max) then the result is undefined.
- If a NaN value is passed to the float functions, NaN is returned.
- By order of performance, for float types, IE is fastest, then EI, then EE, then II. For integer
types, it happens that the order is II, IE and EI ex-aequo, and EE, but the difference is much
less pronounced than for the float types. Use IE preferably for floats, II for ints.
The floating-point functions always compute a random number with the maximum digits of precision,
taking care of bounds. That is, it uses the 1-2 interval as this is a power-of-two bounded interval,
so each representable number in that interval (matching a distinct bit pattern) is given exactly
the same weight. This is really a truly uniform distribution, unlike the 0-1 range (see note below).
That 1-2 interval is then converted to the min-max range, hopefully resulting in the loss of as
few random bits as possible. See also the additional functions below for better performance.
Note: Getting numbers in the 0-1 interval may be tricky. Half the 2^X exponents are in that range as
well as denormal numbers. This could result in an horrible loss of bits when scaling from one exponent
to another. Fortunately, the numbers in 1-2 are generated with maximum precision, and then subtracting 1.0
to get a number in 0-1 keeps that precision because all the numbers obtained this way are still perfectly
representable. Moreover, the minimum exponent reached this way is still far above the denormals range.
Note3: The random number generator MUST be initialized for this function to work correctly.
Note4: These functions are thread-safe if you use one RandomState object per thread (or if you
synchronize the access to a shared state object, of course).
*/
template<bool include_min, bool include_max, typename a_type> a_type Random(a_type min, a_type max, RandomState& state = DefaultRandomState);
// function that returns a random number on the whole possible range
template<typename a_type> a_type Random(RandomState& state = DefaultRandomState);
// Alias that can be useful too
template<typename a_type> inline a_type RandomIE(a_type min, a_type max, RandomState& state = DefaultRandomState) {return Random<true, false, a_type>(min, max, state);}
template<typename a_type> inline a_type RandomEI(a_type min, a_type max, RandomState& state = DefaultRandomState) {return Random<false, true, a_type>(min, max, state);}
template<typename a_type> inline a_type RandomEE(a_type min, a_type max, RandomState& state = DefaultRandomState) {return Random<false, false, a_type>(min, max, state);}
template<typename a_type> inline a_type RandomII(a_type min, a_type max, RandomState& state = DefaultRandomState) {return Random<true, true, a_type>(min, max, state);}
#define STREFLOP_RANDOM_MAKE_REAL(a_type) \
template<> a_type Random<a_type>(RandomState& state); \
template<> a_type Random<true, true, a_type>(a_type min, a_type max, RandomState& state); \
template<> a_type Random<true, false, a_type>(a_type min, a_type max, RandomState& state); \
template<> a_type Random<false, true, a_type>(a_type min, a_type max, RandomState& state); \
template<> a_type Random<false, false, a_type>(a_type min, a_type max, RandomState& state);
STREFLOP_RANDOM_MAKE_REAL(char)
STREFLOP_RANDOM_MAKE_REAL(unsigned char)
STREFLOP_RANDOM_MAKE_REAL(short)
STREFLOP_RANDOM_MAKE_REAL(unsigned short)
STREFLOP_RANDOM_MAKE_REAL(int)
STREFLOP_RANDOM_MAKE_REAL(unsigned int)
STREFLOP_RANDOM_MAKE_REAL(long)
STREFLOP_RANDOM_MAKE_REAL(unsigned long)
STREFLOP_RANDOM_MAKE_REAL(long long)
STREFLOP_RANDOM_MAKE_REAL(unsigned long long)
/** Additional and faster functions for real numbers
These return a number in the 1..2 range, the base for all the other random functions
*/
template<bool include_min, bool include_max, typename a_type> a_type Random12(RandomState& state = DefaultRandomState);
// Alias that can be useful too
template<typename a_type> inline a_type Random12IE(RandomState& state = DefaultRandomState) {return Random12<true, false, a_type>(state);}
template<typename a_type> inline a_type Random12EI(RandomState& state = DefaultRandomState) {return Random12<false, true, a_type>(state);}
template<typename a_type> inline a_type Random12EE(RandomState& state = DefaultRandomState) {return Random12<false, false, a_type>(state);}
template<typename a_type> inline a_type Random12II(RandomState& state = DefaultRandomState) {return Random12<true, true, a_type>(state);}
/** Additional and faster functions for real numbers
These return a number in the 0..1 range by subtracting one to the 1..2 function
This avoids a multiplication for the min...max scaling
Provided for convenience only, use the 1..2 range as the base random function
*/
template<bool include_min, bool include_max, typename a_type> inline a_type Random01(RandomState& state = DefaultRandomState) {
return Random12<include_min, include_max, a_type>(state) - a_type(1.0);
}
// Alias that can be useful too
template<typename a_type> inline a_type Random01IE(RandomState& state = DefaultRandomState) {return Random01<true, false, a_type>(state);}
template<typename a_type> inline a_type Random01EI(RandomState& state = DefaultRandomState) {return Random01<false, true, a_type>(state);}
template<typename a_type> inline a_type Random01EE(RandomState& state = DefaultRandomState) {return Random01<false, false, a_type>(state);}
template<typename a_type> inline a_type Random01II(RandomState& state = DefaultRandomState) {return Random01<true, true, a_type>(state);}
/// Define all 12 and 01 functions only for real types
/// use the 12 function to generate the other
#define STREFLOP_RANDOM_MAKE_REAL_FLOAT_TYPES(a_type) \
template<> a_type Random12<true, true, a_type>(RandomState& state); \
template<> a_type Random12<true, false, a_type>(RandomState& state); \
template<> a_type Random12<false, true, a_type>(RandomState& state); \
template<> a_type Random12<false, false, a_type>(RandomState& state); \
template<> a_type Random<a_type>(RandomState& state); \
template<> inline a_type Random<true, true, a_type>(a_type min, a_type max, RandomState& state) { \
a_type range = max - min;\
return Random12<true,true,a_type>(DefaultRandomState) * range - range + min;\
} \
template<> inline a_type Random<true, false, a_type>(a_type min, a_type max, RandomState& state) { \
a_type range = max - min;\
return Random12<true,false,a_type>(DefaultRandomState) * range - range + min;\
} \
template<> inline a_type Random<false, true, a_type>(a_type min, a_type max, RandomState& state) { \
a_type range = max - min;\
return Random12<false,true,a_type>(DefaultRandomState) * range - range + min;\
} \
template<> inline a_type Random<false, false, a_type>(a_type min, a_type max, RandomState& state) { \
a_type range = max - min;\
return Random12<false,false,a_type>(DefaultRandomState) * range - range + min;\
}
STREFLOP_RANDOM_MAKE_REAL_FLOAT_TYPES(Simple)
STREFLOP_RANDOM_MAKE_REAL_FLOAT_TYPES(Double)
#if defined(Extended)
STREFLOP_RANDOM_MAKE_REAL_FLOAT_TYPES(Extended)
#endif
/**
Utility to get a number from a Normal distribution
The no argument version returns a distribution with mean 0, variance 1.
The 2 argument version takes the desired mean standard deviation.
Both are only defined over the floating-point types
This uses the common polar-coordinate method to transform between uniform and normal
Step 1: generate a pair of numbers in the -1,1 x -1,1 square
Step 2: loop over to step 1 until it is also in the unit circle
Step 3: Convert using the property property (U1, U2) = (X,Y) * sqrt( -2 * log(D) / D)
with D = X*X+Y*Y the squared distance and log(D) the natural logarithm function
Step 4: Choose U1 or U2, they are independent (keep the other for the next call)
Step 5 (optional): Scale and translate U to a given mean, std_dev
There may be better numerical methods, but I'm too lazy to implement them.
Any suggestion/contribution is welcome!
In particular, it seems quite horrendous to do computations close to the 0, in the 0-1 range,
where there may be denormals and the lot, to scale and translate later on. It would be better
to compute everything in another float interval (ex: 2 to 4 centered on 3 for the same square
has a uniform representable range of floats), and only then convert to the given mean/std_dev
at the last moment.
Note: An optional argument "secondary" may be specified, and in that case, a second number
indepedent from the first will be returned at negligible cost (in other words, by default, half the values
are thrown away)
*/
template<typename a_type> a_type NRandom(a_type mean, a_type std_dev, a_type *secondary = 0, RandomState& state = DefaultRandomState);
template<> Simple NRandom(Simple mean, Simple std_dev, Simple *secondary, RandomState& state);
template<> Double NRandom(Double mean, Double std_dev, Double *secondary, RandomState& state);
#if defined(Extended)
template<> Extended NRandom(Extended mean, Extended std_dev, Extended *secondary, RandomState& state);
#endif
/// Simplified versions
template<typename a_type> a_type NRandom(a_type *secondary = 0, RandomState& state = DefaultRandomState);
template<> Simple NRandom(Simple *secondary, RandomState& state);
template<> Double NRandom(Double *secondary, RandomState& state);
#if defined(Extended)
template<> Extended NRandom(Extended *secondary, RandomState& state);
#endif
}
#endif

721
source/shared_lib/include/streflop/SMath.h Normal file
View File

@@ -0,0 +1,721 @@
/*
streflop: STandalone REproducible FLOating-Point
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
Heavily relies on GNU Libm, itself depending on netlib fplibm, GNU MP,
and IBM MP lib.
Uses SoftFloat too.
Please read the history and copyright information in the documentation
provided with the source code
*/
// Included by the main streflop include file
// module broken apart for logical code separation
#ifndef STREFLOP_MATH_H
#define STREFLOP_MATH_H
// just in case, should already be included
#include "streflop.h"
// Names from the libm conversion
namespace streflop_libm {
using streflop::Simple;
using streflop::Double;
#ifdef Extended
using streflop::Extended;
#endif
extern Simple __ieee754_sqrtf(Simple x);
extern Simple __cbrtf(Simple x);
extern Simple __ieee754_hypotf(Simple x, Simple y);
extern Simple __ieee754_expf(Simple x);
extern Simple __ieee754_logf(Simple x);
extern Simple __ieee754_log2f(Simple x);
extern Simple __ieee754_exp2f(Simple x);
extern Simple __ieee754_log10f(Simple x);
extern Simple __ieee754_powf(Simple x, Simple y);
extern Simple __sinf(Simple x);
extern Simple __cosf(Simple x);
extern Simple __tanhf(Simple x);
extern Simple __tanf(Simple x);
extern Simple __ieee754_acosf(Simple x);
extern Simple __ieee754_asinf(Simple x);
extern Simple __atanf(Simple x);
extern Simple __ieee754_atan2f(Simple x, Simple y);
extern Simple __ieee754_coshf(Simple x);
extern Simple __ieee754_sinhf(Simple x);
extern Simple __ieee754_acoshf(Simple x);
extern Simple __asinhf(Simple x);
extern Simple __ieee754_atanhf(Simple x);
extern Simple __fabsf(Simple x);
extern Simple __floorf(Simple x);
extern Simple __ceilf(Simple x);
extern Simple __truncf(Simple x);
extern Simple __ieee754_fmodf(Simple x, Simple y);
extern Simple __ieee754_remainderf(Simple x, Simple y);
extern Simple __remquof(Simple x, Simple y, int *quo);
extern Simple __rintf(Simple x);
extern long int __lrintf(Simple x);
extern long long int __llrintf(Simple x);
extern Simple __roundf(Simple x);
extern long int __lroundf(Simple x);
extern long long int __llroundf(Simple x);
extern Simple __nearbyintf(Simple x);
extern Simple __frexpf(Simple x, int *exp);
extern Simple __ldexpf(Simple value, int exp);
extern Simple __logbf(Simple x);
extern int __ilogbf(Simple x);
extern Simple __copysignf(Simple x);
extern int __signbitf(Simple x);
extern Simple __nextafterf(Simple x, Simple y);
extern Simple __expm1f(Simple x);
extern Simple __log1pf(Simple x);
extern Simple __erff(Simple x);
extern Simple __ieee754_j0f(Simple x);
extern Simple __ieee754_j1f(Simple x);
extern Simple __ieee754_jnf(int n, Simple x);
extern Simple __ieee754_y0f(Simple x);
extern Simple __ieee754_y1f(Simple x);
extern Simple __ieee754_ynf(int n, Simple x);
extern Simple __scalbnf(Simple x, int n);
extern Simple __scalblnf(Simple x, long int n);
extern int __fpclassifyf(Simple x);
extern int __isnanf(Simple x);
extern int __isinff(Simple x);
extern Double __ieee754_sqrt(Double x);
extern Double __cbrt(Double x);
extern Double __ieee754_hypot(Double x, Double y);
extern Double __ieee754_exp(Double x);
extern Double __ieee754_log(Double x);
extern Double __ieee754_log2(Double x);
extern Double __ieee754_exp2(Double x);
extern Double __ieee754_log10(Double x);
extern Double __ieee754_pow(Double x, Double y);
extern Double __sin(Double x);
extern Double __cos(Double x);
extern Double tan(Double x);
extern Double __ieee754_acos(Double x);
extern Double __ieee754_asin(Double x);
extern Double atan(Double x);
extern Double __ieee754_atan2(Double x, Double y);
extern Double __ieee754_cosh(Double x);
extern Double __ieee754_sinh(Double x);
extern Double __tanh(Double x);
extern Double __ieee754_acosh(Double x);
extern Double __asinh(Double x);
extern Double __ieee754_atanh(Double x);
extern Double __fabs(Double x);
extern Double __floor(Double x);
extern Double __ceil(Double x);
extern Double __trunc(Double x);
extern Double __ieee754_fmod(Double x, Double y);
extern Double __ieee754_remainder(Double x, Double y);
extern Double __remquo(Double x, Double y, int *quo);
extern Double __rint(Double x);
extern long int __lrint(Double x);
extern long long int __llrint(Double x);
extern Double __round(Double x);
extern long int __lround(Double x);
extern long long int __llround(Double x);
extern Double __nearbyint(Double x);
extern Double __frexp(Double x, int *exp);
extern Double __ldexp(Double value, int exp);
extern Double __logb(Double x);
extern int __ilogb(Double x);
extern Double __copysign(Double x);
extern int __signbit(Double x);
extern Double __nextafter(Double x, Double y);
extern Double __expm1(Double x);
extern Double __log1p(Double x);
extern Double __erf(Double x);
extern Double __ieee754_j0(Double x);
extern Double __ieee754_j1(Double x);
extern Double __ieee754_jn(int n, Double x);
extern Double __ieee754_y0(Double x);
extern Double __ieee754_y1(Double x);
extern Double __ieee754_yn(int n, Double x);
extern Double __scalbn(Double x, int n);
extern Double __scalbln(Double x, long int n);
extern int __fpclassify(Double x);
extern int __isnanl(Double x);
extern int __isinf(Double x);
#ifdef Extended
extern Extended __ieee754_sqrtl(Extended x);
extern Extended __cbrtl(Extended x);
extern Extended __ieee754_hypotl(Extended x, Extended y);
extern Extended __ieee754_expl(Extended x);
extern Extended __ieee754_logl(Extended x);
extern Extended __sinl(Extended x);
extern Extended __cosl(Extended x);
extern Extended __tanl(Extended x);
extern Extended __ieee754_asinl(Extended x);
extern Extended __atanl(Extended x);
extern Extended __ieee754_atan2l(Extended x, Extended y);
extern Extended __ieee754_coshl(Extended x);
extern Extended __ieee754_sinhl(Extended x);
extern Extended __tanhl(Extended x);
extern Extended __ieee754_acoshl(Extended x);
extern Extended __asinhl(Extended x);
extern Extended __ieee754_atanhl(Extended x);
extern Extended __fabsl(Extended x);
extern Extended __floorl(Extended x);
extern Extended __ceill(Extended x);
extern Extended __truncl(Extended x);
extern Extended __ieee754_fmodl(Extended x, Extended y);
extern Extended __ieee754_remainderl(Extended x, Extended y);
extern Extended __remquol(Extended x, Extended y, int *quo);
extern Extended __rintl(Extended x);
extern long int __lrintl(Extended x);
extern long long int __llrintl(Extended x);
extern Extended __roundl(Extended x);
extern long int __lroundl(Extended x);
extern long long int __llroundl(Extended x);
extern Extended __nearbyintl(Extended x);
extern Extended __frexpl(Extended x, int *exp);
extern Extended __ldexpl(Extended value, int exp);
extern Extended __logbl(Extended x);
extern int __ilogbl(Extended x);
extern Extended __copysignl(Extended x);
extern int __signbitl(Extended x);
extern Extended __nextafterl(Extended x, Extended y);
extern Extended __expm1l(Extended x);
extern Extended __log1pl(Extended x);
extern Extended __erfl(Extended x);
extern Extended __ieee754_j0l(Extended x);
extern Extended __ieee754_j1l(Extended x);
extern Extended __ieee754_jnl(int n, Extended x);
extern Extended __ieee754_y0l(Extended x);
extern Extended __ieee754_y1l(Extended x);
extern Extended __ieee754_ynl(int n, Extended x);
extern Extended __scalbnl(Extended x, int n);
extern Extended __scalblnl(Extended x, long int n);
extern int __fpclassifyl(Extended x);
extern int __isnanl(Extended x);
extern int __isinfl(Extended x);
#endif // Extended
}
// Wrappers in our own namespace
namespace streflop {
// Stolen from math.h. All floating-point numbers can be put in one of these categories.
enum
{
STREFLOP_FP_NAN = 0,
STREFLOP_FP_INFINITE = 1,
STREFLOP_FP_ZERO = 2,
STREFLOP_FP_SUBNORMAL = 3,
STREFLOP_FP_NORMAL = 4
};
// Simple and double are present in all configurations
inline Simple sqrt(Simple x) {return streflop_libm::__ieee754_sqrtf(x);}
inline Simple cbrt(Simple x) {return streflop_libm::__cbrtf(x);}
inline Simple hypot(Simple x, Simple y) {return streflop_libm::__ieee754_hypotf(x,y);}
inline Simple exp(Simple x) {return streflop_libm::__ieee754_expf(x);}
inline Simple log(Simple x) {return streflop_libm::__ieee754_logf(x);}
inline Simple log2(Simple x) {return streflop_libm::__ieee754_log2f(x);}
inline Simple exp2(Simple x) {return streflop_libm::__ieee754_exp2f(x);}
inline Simple log10(Simple x) {return streflop_libm::__ieee754_log10f(x);}
inline Simple pow(Simple x, Simple y) {return streflop_libm::__ieee754_powf(x,y);}
inline Simple sin(Simple x) {return streflop_libm::__sinf(x);}
inline Simple cos(Simple x) {return streflop_libm::__cosf(x);}
inline Simple tan(Simple x) {return streflop_libm::__tanf(x);}
inline Simple acos(Simple x) {return streflop_libm::__ieee754_acosf(x);}
inline Simple asin(Simple x) {return streflop_libm::__ieee754_asinf(x);}
inline Simple atan(Simple x) {return streflop_libm::__atanf(x);}
inline Simple atan2(Simple x, Simple y) {return streflop_libm::__ieee754_atan2f(x,y);}
inline Simple cosh(Simple x) {return streflop_libm::__ieee754_coshf(x);}
inline Simple sinh(Simple x) {return streflop_libm::__ieee754_sinhf(x);}
inline Simple tanh(Simple x) {return streflop_libm::__tanhf(x);}
inline Simple acosh(Simple x) {return streflop_libm::__ieee754_acoshf(x);}
inline Simple asinh(Simple x) {return streflop_libm::__asinhf(x);}
inline Simple atanh(Simple x) {return streflop_libm::__ieee754_atanhf(x);}
inline Simple fabs(Simple x) {return streflop_libm::__fabsf(x);}
inline Simple floor(Simple x) {return streflop_libm::__floorf(x);}
inline Simple ceil(Simple x) {return streflop_libm::__ceilf(x);}
inline Simple trunc(Simple x) {return streflop_libm::__truncf(x);}
inline Simple fmod(Simple x, Simple y) {return streflop_libm::__ieee754_fmodf(x,y);}
inline Simple remainder(Simple x, Simple y) {return streflop_libm::__ieee754_remainderf(x,y);}
inline Simple remquo(Simple x, Simple y, int *quo) {return streflop_libm::__remquof(x,y,quo);}
inline Simple rint(Simple x) {return streflop_libm::__rintf(x);}
inline long int lrint(Simple x) {return streflop_libm::__lrintf(x);}
inline long long int llrint(Simple x) {return streflop_libm::__llrintf(x);}
inline Simple round(Simple x) {return streflop_libm::__roundf(x);}
inline long int lround(Simple x) {return streflop_libm::__lroundf(x);}
inline long long int llround(Simple x) {return streflop_libm::__llroundf(x);}
inline Simple nearbyint(Simple x) {return streflop_libm::__nearbyintf(x);}
inline Simple frexp(Simple x, int *exp) {return streflop_libm::__frexpf(x,exp);}
inline Simple ldexp(Simple value, int exp) {return streflop_libm::__ldexpf(value,exp);}
inline Simple logb(Simple x) {return streflop_libm::__logbf(x);}
inline int ilogb(Simple x) {return streflop_libm::__ilogbf(x);}
inline Simple copysign(Simple x) {return streflop_libm::__copysignf(x);}
#undef signbit
inline int signbit (Simple x) {return streflop_libm::__signbitf(x);}
inline Simple nextafter(Simple x, Simple y) {return streflop_libm::__nextafterf(x,y);}
inline Simple expm1(Simple x) {return streflop_libm::__expm1f(x);}
inline Simple log1p(Simple x) {return streflop_libm::__log1pf(x);}
inline Simple erf(Simple x) {return streflop_libm::__erff(x);}
inline Simple j0(Simple x) {return streflop_libm::__ieee754_j0f(x);}
inline Simple j1(Simple x) {return streflop_libm::__ieee754_j1f(x);}
inline Simple jn(int n, Simple x) {return streflop_libm::__ieee754_jnf(n,x);}
inline Simple y0(Simple x) {return streflop_libm::__ieee754_y0f(x);}
inline Simple y1(Simple x) {return streflop_libm::__ieee754_y1f(x);}
inline Simple yn(int n, Simple x) {return streflop_libm::__ieee754_ynf(n,x);}
inline Simple scalbn(Simple x, int n) {return streflop_libm::__scalbnf(x,n);}
inline Simple scalbln(Simple x, long int n) {return streflop_libm::__scalblnf(x,n);}
#undef fpclassify
inline int fpclassify(Simple x) {return streflop_libm::__fpclassifyf(x);}
#undef isnan
inline int isnan(Simple x) {return streflop_libm::__isnanf(x);}
#undef isinf
inline int isinf(Simple x) {return streflop_libm::__isinff(x);}
#undef isfinite
inline int isfinite(Simple x) {return !(isnan(x) || isinf(x));}
// Stolen from math.h and inlined instead of macroized.
// Return nonzero value if X is neither zero, subnormal, Inf, nor NaN. */
#undef isnormal
inline int isnormal(Simple x) {return fpclassify(x) == STREFLOP_FP_NORMAL;}
// Constants user may set
extern const Simple SimplePositiveInfinity;
extern const Simple SimpleNegativeInfinity;
// Non-signaling NaN used for returning such results
// Standard lets a large room for implementing different kinds of NaN
// These NaN can be used for custom purposes
// Threated as an integer, the bit pattern here may be incremented to get at least 2^20 different NaN custom numbers!
// Note that when switching the left-most bit to 1, you can get another bunch of negative NaNs, whatever this mean.
extern const Simple SimpleNaN;
/** Generic C99 "macros" for unordered comparison
Defined as inlined for each type, thanks to C++ overloading
*/
#if defined isunordered
#undef isunordered
#endif // defined isunordered
inline bool isunordered(Simple x, Simple y) {
return (fpclassify(x) == STREFLOP_FP_NAN) || (fpclassify(y) == STREFLOP_FP_NAN);
}
#if defined isgreater
#undef isgreater
#endif // defined isgreater
inline bool isgreater(Simple x, Simple y) {
return (!isunordered(x,y)) && (x > y);
}
#if defined isgreaterequal
#undef isgreaterequal
#endif // defined isgreaterequal
inline bool isgreaterequal(Simple x, Simple y) {
return (!isunordered(x,y)) && (x >= y);
}
#if defined isless
#undef isless
#endif // defined isless
inline bool isless(Simple x, Simple y) {
return (!isunordered(x,y)) && (x < y);
}
#if defined islessequal
#undef islessequal
#endif // defined islessequal
inline bool islessequal(Simple x, Simple y) {
return (!isunordered(x,y)) && (x <= y);
}
#if defined islessgreater
#undef islessgreater
#endif // defined islessgreater
inline bool islessgreater(Simple x, Simple y) {
return (!isunordered(x,y)) && ((x < y) || (x > y));
}
// Add xxxf alias to ease porting existing code to streflop
// Additionally, using xxxf(number) avoids potential confusion
inline Simple sqrtf(Simple x) {return sqrt(x);}
inline Simple cbrtf(Simple x) {return cbrt(x);}
inline Simple hypotf(Simple x, Simple y) {return hypot(x, y);}
inline Simple expf(Simple x) {return exp(x);}
inline Simple logf(Simple x) {return log(x);}
inline Simple log2f(Simple x) {return log2(x);}
inline Simple exp2f(Simple x) {return exp2(x);}
inline Simple log10f(Simple x) {return log10(x);}
inline Simple powf(Simple x, Simple y) {return pow(x, y);}
inline Simple sinf(Simple x) {return sin(x);}
inline Simple cosf(Simple x) {return cos(x);}
inline Simple tanf(Simple x) {return tan(x);}
inline Simple acosf(Simple x) {return acos(x);}
inline Simple asinf(Simple x) {return asin(x);}
inline Simple atanf(Simple x) {return atan(x);}
inline Simple atan2f(Simple x, Simple y) {return atan2(x, y);}
inline Simple coshf(Simple x) {return cosh(x);}
inline Simple sinhf(Simple x) {return sinh(x);}
inline Simple tanhf(Simple x) {return tanh(x);}
inline Simple acoshf(Simple x) {return acosh(x);}
inline Simple asinhf(Simple x) {return asinh(x);}
inline Simple atanhf(Simple x) {return atanh(x);}
inline Simple fabsf(Simple x) {return fabs(x);}
inline Simple floorf(Simple x) {return floor(x);}
inline Simple ceilf(Simple x) {return ceil(x);}
inline Simple truncf(Simple x) {return trunc(x);}
inline Simple fmodf(Simple x, Simple y) {return fmod(x,y);}
inline Simple remainderf(Simple x, Simple y) {return remainder(x,y);}
inline Simple remquof(Simple x, Simple y, int *quo) {return remquo(x, y, quo);}
inline Simple rintf(Simple x) {return rint(x);}
inline long int lrintf(Simple x) {return lrint(x);}
inline long long int llrintf(Simple x) {return llrint(x);}
inline Simple roundf(Simple x) {return round(x);}
inline long int lroundf(Simple x) {return lround(x);}
inline long long int llroundf(Simple x) {return llround(x);}
inline Simple nearbyintf(Simple x) {return nearbyint(x);}
inline Simple frexpf(Simple x, int *exp) {return frexp(x, exp);}
inline Simple ldexpf(Simple value, int exp) {return ldexp(value,exp);}
inline Simple logbf(Simple x) {return logb(x);}
inline int ilogbf(Simple x) {return ilogb(x);}
inline Simple copysignf(Simple x) {return copysign(x);}
inline int signbitf(Simple x) {return signbit(x);}
inline Simple nextafterf(Simple x, Simple y) {return nextafter(x, y);}
inline Simple expm1f(Simple x) {return expm1(x);}
inline Simple log1pf(Simple x) {return log1p(x);}
inline Simple erff(Simple x) {return erf(x);}
inline Simple j0f(Simple x) {return j0(x);}
inline Simple j1f(Simple x) {return j1(x);}
inline Simple jnf(int n, Simple x) {return jn(n, x);}
inline Simple y0f(Simple x) {return y0(x);}
inline Simple y1f(Simple x) {return y1(x);}
inline Simple ynf(int n, Simple x) {return yn(n, x);}
inline Simple scalbnf(Simple x, int n) {return scalbn(x, n);}
inline Simple scalblnf(Simple x, long int n) {return scalbln(x, n);}
inline int fpclassifyf(Simple x) {return fpclassify(x);}
inline int isnanf(Simple x) {return isnan(x);}
inline int isinff(Simple x) {return isinf(x);}
inline int isfinitef(Simple x) {return isfinite(x);}
inline int isnormalf(Simple x) {return isnormal(x);}
inline bool isunorderedf(Simple x, Simple y) {return isunordered(x, y);}
inline bool isgreaterf(Simple x, Simple y) {return isgreater(x, y);}
inline bool isgreaterequalf(Simple x, Simple y) {return isgreaterequalf(x, y);}
inline bool islessf(Simple x, Simple y) {return isless(x, y);}
inline bool islessequalf(Simple x, Simple y) {return islessequal(x, y);}
inline bool islessgreaterf(Simple x, Simple y) {return islessgreater(x, y);}
// Declare Double functions
// Simple and double are present in all configurations
inline Double sqrt(Double x) {return streflop_libm::__ieee754_sqrt(x);}
inline Double cbrt(Double x) {return streflop_libm::__cbrt(x);}
inline Double hypot(Double x, Double y) {return streflop_libm::__ieee754_hypot(x,y);}
inline Double exp(Double x) {return streflop_libm::__ieee754_exp(x);}
inline Double log(Double x) {return streflop_libm::__ieee754_log(x);}
inline Double log2(Double x) {return streflop_libm::__ieee754_log2(x);}
inline Double exp2(Double x) {return streflop_libm::__ieee754_exp2(x);}
inline Double log10(Double x) {return streflop_libm::__ieee754_log10(x);}
inline Double pow(Double x, Double y) {return streflop_libm::__ieee754_pow(x,y);}
inline Double sin(Double x) {return streflop_libm::__sin(x);}
inline Double cos(Double x) {return streflop_libm::__cos(x);}
inline Double tan(Double x) {return streflop_libm::tan(x);}
inline Double acos(Double x) {return streflop_libm::__ieee754_acos(x);}
inline Double asin(Double x) {return streflop_libm::__ieee754_asin(x);}
inline Double atan(Double x) {return streflop_libm::atan(x);}
inline Double atan2(Double x, Double y) {return streflop_libm::__ieee754_atan2(x,y);}
inline Double cosh(Double x) {return streflop_libm::__ieee754_cosh(x);}
inline Double sinh(Double x) {return streflop_libm::__ieee754_sinh(x);}
inline Double tanh(Double x) {return streflop_libm::__tanh(x);}
inline Double acosh(Double x) {return streflop_libm::__ieee754_acosh(x);}
inline Double asinh(Double x) {return streflop_libm::__asinh(x);}
inline Double atanh(Double x) {return streflop_libm::__ieee754_atanh(x);}
inline Double fabs(Double x) {return streflop_libm::__fabs(x);}
inline Double floor(Double x) {return streflop_libm::__floor(x);}
inline Double ceil(Double x) {return streflop_libm::__ceil(x);}
inline Double trunc(Double x) {return streflop_libm::__trunc(x);}
inline Double fmod(Double x, Double y) {return streflop_libm::__ieee754_fmod(x,y);}
inline Double remainder(Double x, Double y) {return streflop_libm::__ieee754_remainder(x,y);}
inline Double remquo(Double x, Double y, int *quo) {return streflop_libm::__remquo(x,y,quo);}
inline Double rint(Double x) {return streflop_libm::__rint(x);}
inline long int lrint(Double x) {return streflop_libm::__lrint(x);}
inline long long int llrint(Double x) {return streflop_libm::__llrint(x);}
inline Double round(Double x) {return streflop_libm::__round(x);}
inline long int lround(Double x) {return streflop_libm::__lround(x);}
inline long long int llround(Double x) {return streflop_libm::__llround(x);}
inline Double nearbyint(Double x) {return streflop_libm::__nearbyint(x);}
inline Double frexp(Double x, int *exp) {return streflop_libm::__frexp(x, exp);}
inline Double ldexp(Double value, int exp) {return streflop_libm::__ldexp(value,exp);}
inline Double logb(Double x) {return streflop_libm::__logb(x);}
inline int ilogb(Double x) {return streflop_libm::__ilogb(x);}
inline Double copysign(Double x) {return streflop_libm::__copysign(x);}
inline int signbit(Double x) {return streflop_libm::__signbit(x);}
inline Double nextafter(Double x, Double y) {return streflop_libm::__nextafter(x,y);}
inline Double expm1(Double x) {return streflop_libm::__expm1(x);}
inline Double log1p(Double x) {return streflop_libm::__log1p(x);}
inline Double erf(Double x) {return streflop_libm::__erf(x);}
inline Double j0(Double x) {return streflop_libm::__ieee754_j0(x);}
inline Double j1(Double x) {return streflop_libm::__ieee754_j1(x);}
inline Double jn(int n, Double x) {return streflop_libm::__ieee754_jn(n,x);}
inline Double y0(Double x) {return streflop_libm::__ieee754_y0(x);}
inline Double y1(Double x) {return streflop_libm::__ieee754_y1(x);}
inline Double yn(int n, Double x) {return streflop_libm::__ieee754_yn(n,x);}
inline Double scalbn(Double x, int n) {return streflop_libm::__scalbn(x,n);}
inline Double scalbln(Double x, long int n) {return streflop_libm::__scalbln(x,n);}
inline int fpclassify(Double x) {return streflop_libm::__fpclassify(x);}
inline int isnan(Double x) {return streflop_libm::__isnanl(x);}
inline int isinf(Double x) {return streflop_libm::__isinf(x);}
inline int isfinite(Double x) {return !(isnan(x) || isinf(x));}
// Stolen from math.h and inlined instead of macroized.
// Return nonzero value if X is neither zero, subnormal, Inf, nor NaN. */
inline int isnormal(Double x) {return fpclassify(x) == STREFLOP_FP_NORMAL;}
// Constants user may set
extern const Double DoublePositiveInfinity;
extern const Double DoubleNegativeInfinity;
// Non-signaling NaN used for returning such results
// Standard lets a large room for implementing different kinds of NaN
// These NaN can be used for custom purposes
// Threated as an integer, the bit pattern here may be incremented to get at least 2^20 different NaN custom numbers!
// Note that when switching the left-most bit to 1, you can get another bunch of negative NaNs, whatever this mean.
extern const Double DoubleNaN;
/** Generic C99 "macros" for unordered comparison
Defined as inlined for each type, thanks to C++ overloading
*/
inline bool isunordered(Double x, Double y) {
return (fpclassify(x) == STREFLOP_FP_NAN) || (fpclassify (y) == STREFLOP_FP_NAN);
}
inline bool isgreater(Double x, Double y) {
return (!isunordered(x,y)) && (x > y);
}
inline bool isgreaterequal(Double x, Double y) {
return (!isunordered(x,y)) && (x >= y);
}
inline bool isless(Double x, Double y) {
return (!isunordered(x,y)) && (x < y);
}
inline bool islessequal(Double x, Double y) {
return (!isunordered(x,y)) && (x <= y);
}
inline bool islessgreater(Double x, Double y) {
return (!isunordered(x,y)) && ((x < y) || (x > y));
}
// Extended are not always available
#ifdef Extended
inline Extended cbrt(Extended x) {return streflop_libm::__cbrtl(x);}
inline Extended hypot(Extended x, Extended y) {return streflop_libm::__ieee754_hypotl(x,y);}
// Missing from libm: temporarily switch to Double and execute the Double version,
// then switch back to Extended and return the result
inline Extended sqrt(Extended x) {streflop_init<Double>(); Double res = sqrt(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended exp(Extended x) {streflop_init<Double>(); Double res = exp(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended log(Extended x) {streflop_init<Double>(); Double res = log(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended log2(Extended x) {streflop_init<Double>(); Double res = log2(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended exp2(Extended x) {streflop_init<Double>(); Double res = exp2(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended log10(Extended x) {streflop_init<Double>(); Double res = log10(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended pow(Extended x, Extended y) {streflop_init<Double>(); Double res = pow(Double(x), Double(y)); streflop_init<Extended>(); return Extended(res);}
inline Extended sin(Extended x) {streflop_init<Double>(); Double res = sin(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended cos(Extended x) {streflop_init<Double>(); Double res = cos(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended tan(Extended x) {streflop_init<Double>(); Double res = tan(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended acos(Extended x) {streflop_init<Double>(); Double res = acos(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended asin(Extended x) {streflop_init<Double>(); Double res = asin(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended atan(Extended x) {streflop_init<Double>(); Double res = atan(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended atan2(Extended x, Extended y) {streflop_init<Double>(); Double res = atan2(Double(x), Double(y)); streflop_init<Extended>(); return Extended(res);}
inline Extended cosh(Extended x) {streflop_init<Double>(); Double res = cosh(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended sinh(Extended x) {streflop_init<Double>(); Double res = sinh(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended tanh(Extended x) {streflop_init<Double>(); Double res = tanh(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended acosh(Extended x) {streflop_init<Double>(); Double res = acosh(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended asinh(Extended x) {streflop_init<Double>(); Double res = asinh(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended atanh(Extended x) {streflop_init<Double>(); Double res = atanh(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended expm1(Extended x) {streflop_init<Double>(); Double res = expm1(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended log1p(Extended x) {streflop_init<Double>(); Double res = log1p(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended erf(Extended x) {streflop_init<Double>(); Double res = erf(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended j0(Extended x) {streflop_init<Double>(); Double res = j0(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended j1(Extended x) {streflop_init<Double>(); Double res = j1(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended jn(int n, Extended x) {streflop_init<Double>(); Double res = jn(n,Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended y0(Extended x) {streflop_init<Double>(); Double res = y0(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended y1(Extended x) {streflop_init<Double>(); Double res = y1(Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended yn(int n, Extended x) {streflop_init<Double>(); Double res = yn(n,Double(x)); streflop_init<Extended>(); return Extended(res);}
inline Extended scalbn(Extended x, int n) {streflop_init<Double>(); Double res = scalbn(Double(x),n); streflop_init<Extended>(); return Extended(res);}
inline Extended scalbln(Extended x, long int n) {streflop_init<Double>(); Double res = scalbln(Double(x),n); streflop_init<Extended>(); return Extended(res);}
inline Extended fabs(Extended x) {return streflop_libm::__fabsl(x);}
inline Extended floor(Extended x) {return streflop_libm::__floorl(x);}
inline Extended ceil(Extended x) {return streflop_libm::__ceill(x);}
inline Extended trunc(Extended x) {return streflop_libm::__truncl(x);}
inline Extended fmod(Extended x, Extended y) {return streflop_libm::__ieee754_fmodl(x,y);}
inline Extended remainder(Extended x, Extended y) {return streflop_libm::__ieee754_remainderl(x,y);}
inline Extended remquo(Extended x, Extended y, int *quo) {return streflop_libm::__remquol(x,y,quo);}
inline Extended rint(Extended x) {return streflop_libm::__rintl(x);}
inline long int lrint(Extended x) {return streflop_libm::__lrintl(x);}
inline long long int llrint(Extended x) {return streflop_libm::__llrintl(x);}
inline Extended round(Extended x) {return streflop_libm::__roundl(x);}
inline long int lround(Extended x) {return streflop_libm::__lroundl(x);}
inline long long int llround(Extended x) {return streflop_libm::__llroundl(x);}
inline Extended nearbyint(Extended x) {return streflop_libm::__nearbyintl(x);}
inline Extended frexp(Extended x, int *exp) {return streflop_libm::__frexpl(x,exp);}
inline Extended ldexp(Extended value, int exp) {return streflop_libm::__ldexpl(value,exp);}
inline Extended logb(Extended x) {return streflop_libm::__logbl(x);}
inline int ilogb(Extended x) {return streflop_libm::__ilogbl(x);}
inline Extended copysign(Extended x) {return streflop_libm::__copysignl(x);}
inline int signbit (Extended x) {return streflop_libm::__signbitl(x);}
inline Extended nextafter(Extended x, Extended y) {return streflop_libm::__nextafterl(x,y);}
inline int fpclassify(Extended x) {return streflop_libm::__fpclassifyl(x);}
inline int isnan(Extended x) {return streflop_libm::__isnanl(x);}
inline int isinf(Extended x) {return streflop_libm::__isinfl(x);}
inline int isfinite(Extended x) {return !(isnan(x) || isinf(x));}
// Stolen from math.h and inlined instead of macroized.
// Return nonzero value if X is neither zero, subnormal, Inf, nor NaN. */
inline int isnormal(Extended x) {return fpclassify(x) == STREFLOP_FP_NORMAL;}
// Constants user may set
extern const Extended ExtendedPositiveInfinity;
extern const Extended ExtendedNegativeInfinity;
// Non-signaling NaN used for returning such results
// Standard lets a large room for implementing different kinds of NaN
// These NaN can be used for custom purposes
// Threated as an integer, the bit pattern here may be incremented to get at least 2^20 different NaN custom numbers!
// Note that when switching the left-most bit to 1, you can get another bunch of negative NaNs, whatever this mean.
extern const Extended ExtendedNaN;
/** Generic C99 "macros" for unordered comparison
Defined as inlined for each type, thanks to C++ overloading
*/
inline bool isunordered(Extended x, Extended y) {
return (fpclassify(x) == STREFLOP_FP_NAN) || (fpclassify (y) == STREFLOP_FP_NAN);
}
inline bool isgreater(Extended x, Extended y) {
return (!isunordered(x,y)) && (x > y);
}
inline bool isgreaterequal(Extended x, Extended y) {
return (!isunordered(x,y)) && (x >= y);
}
inline bool isless(Extended x, Extended y) {
return (!isunordered(x,y)) && (x < y);
}
inline bool islessequal(Extended x, Extended y) {
return (!isunordered(x,y)) && (x <= y);
}
inline bool islessgreater(Extended x, Extended y) {
return (!isunordered(x,y)) && ((x < y) || (x > y));
}
// Add xxxl alias to ease porting existing code to streflop
// Additionally, using xxxl(number) avoids potential confusion
inline Extended sqrtl(Extended x) {return sqrt(x);}
inline Extended cbrtl(Extended x) {return cbrt(x);}
inline Extended hypotl(Extended x, Extended y) {return hypot(x, y);}
inline Extended expl(Extended x) {return exp(x);}
inline Extended logl(Extended x) {return log(x);}
inline Extended log2l(Extended x) {return log2(x);}
inline Extended exp2l(Extended x) {return exp2(x);}
inline Extended log10l(Extended x) {return log10(x);}
inline Extended powl(Extended x, Extended y) {return pow(x, y);}
inline Extended sinl(Extended x) {return sin(x);}
inline Extended cosl(Extended x) {return cos(x);}
inline Extended tanl(Extended x) {return tan(x);}
inline Extended acosl(Extended x) {return acos(x);}
inline Extended asinl(Extended x) {return asin(x);}
inline Extended atanl(Extended x) {return atan(x);}
inline Extended atan2l(Extended x, Extended y) {return atan2(x, y);}
inline Extended coshl(Extended x) {return cosh(x);}
inline Extended sinhl(Extended x) {return sinh(x);}
inline Extended tanhl(Extended x) {return tanh(x);}
inline Extended acoshl(Extended x) {return acosh(x);}
inline Extended asinhl(Extended x) {return asinh(x);}
inline Extended atanhl(Extended x) {return atanh(x);}
inline Extended fabsl(Extended x) {return fabs(x);}
inline Extended floorl(Extended x) {return floor(x);}
inline Extended ceill(Extended x) {return ceil(x);}
inline Extended truncl(Extended x) {return trunc(x);}
inline Extended fmodl(Extended x, Extended y) {return fmod(x,y);}
inline Extended remainderl(Extended x, Extended y) {return remainder(x,y);}
inline Extended remquol(Extended x, Extended y, int *quo) {return remquo(x, y, quo);}
inline Extended rintl(Extended x) {return rint(x);}
inline long int lrintl(Extended x) {return lrint(x);}
inline long long int llrintl(Extended x) {return llrint(x);}
inline Extended roundl(Extended x) {return round(x);}
inline long int lroundl(Extended x) {return lround(x);}
inline long long int llroundl(Extended x) {return llround(x);}
inline Extended nearbyintl(Extended x) {return nearbyint(x);}
inline Extended frexpl(Extended x, int *exp) {return frexp(x, exp);}
inline Extended ldexpl(Extended value, int exp) {return ldexp(value,exp);}
inline Extended logbl(Extended x) {return logb(x);}
inline int ilogbl(Extended x) {return ilogb(x);}
inline Extended copysignl(Extended x) {return copysign(x);}
inline int signbitl(Extended x) {return signbit(x);}
inline Extended nextafterl(Extended x, Extended y) {return nextafter(x, y);}
inline Extended expm1l(Extended x) {return expm1(x);}
inline Extended log1pl(Extended x) {return log1p(x);}
inline Extended erfl(Extended x) {return erf(x);}
inline Extended j0l(Extended x) {return j0(x);}
inline Extended j1l(Extended x) {return j1(x);}
inline Extended jnl(int n, Extended x) {return jn(n, x);}
inline Extended y0l(Extended x) {return y0(x);}
inline Extended y1l(Extended x) {return y1(x);}
inline Extended ynl(int n, Extended x) {return yn(n, x);}
inline Extended scalbnl(Extended x, int n) {return scalbn(x, n);}
inline Extended scalblnl(Extended x, long int n) {return scalbln(x, n);}
inline int fpclassifyl(Extended x) {return fpclassify(x);}
inline int isnanl(Extended x) {return isnan(x);}
inline int isinfl(Extended x) {return isinf(x);}
inline int isfinitel(Extended x) {return isfinite(x);}
inline int isnormall(Extended x) {return isnormal(x);}
inline bool isunorderedl(Extended x, Extended y) {return isunordered(x, y);}
inline bool isgreaterl(Extended x, Extended y) {return isgreater(x, y);}
inline bool isgreaterequall(Extended x, Extended y) {return isgreaterequall(x, y);}
inline bool islessl(Extended x, Extended y) {return isless(x, y);}
inline bool islessequall(Extended x, Extended y) {return islessequal(x, y);}
inline bool islessgreaterl(Extended x, Extended y) {return islessgreater(x, y);}
#endif // Extended
} // namespace streflop
#endif // STREFLOP_MATH_H

82
source/shared_lib/include/streflop/System.h Normal file
View File

@@ -0,0 +1,82 @@
/*
streflop: STandalone REproducible FLOating-Point
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
Heavily relies on GNU Libm, itself depending on netlib fplibm, GNU MP, and IBM MP lib.
Uses SoftFloat too.
Please read the history and copyright information in the documentation provided with the source code
*/
#ifndef STREFLOP_SYSTEM_H
#define STREFLOP_SYSTEM_H
#if defined(STREFLOP_SSE) || defined(STREFLOP_X87)
// SSE or X87 machines are little-endian
#define __BYTE_ORDER 1234
#define __FLOAT_WORD_ORDER 1234
// Softfloat or other unknown FPU. TODO: Try some header autodetect?
#else
#define __BYTE_ORDER 1234
#define __FLOAT_WORD_ORDER 1234
/*
#include <endian.h>
#if !defined(__FLOAT_WORD_ORDER) || !defined(__BYTE_ORDER)
#undef __FLOAT_WORD_ORDER
#undef __BYTE_ORDER
#if defined(_LITTLE_ENDIAN) || defined(__LITTLE_ENDIAN) || defined(LITTLE_ENDIAN)
// avoid conflicting defs
#undef _LITTLE_ENDIAN
#undef LITTLE_ENDIAN
#undef __LITTLE_ENDIAN
#undef __BYTE_ORDER
#define __BYTE_ORDER 1234
#define __FLOAT_WORD_ORDER 1234
#elif defined(_BIG_ENDIAN) || defined(__BIG_ENDIAN) || defined(BIG_ENDIAN)
// avoid conflicting defs
#undef _BIG_ENDIAN
#undef BIG_ENDIAN
#undef __BIG_ENDIAN
#undef __BYTE_ORDER
#define __BYTE_ORDER 4321
#define __FLOAT_WORD_ORDER 4321
#else
// Unfortunately, the endianity for the compilation target cannot be deduced at compile time
// => preprocessor #define is necessary, possibly with heuristics
#error Could not detect endianity.
#endif
#endif
*/
#endif
#include "IntegerTypes.h"
// Avoid conflict with system types, if any
namespace streflop {
// Use the assumption above for char
typedef char int8_t;
typedef unsigned char uint8_t;
// Now "run" the meta-program to define the types
// Entry point is (unsigned) short, then all types above are tried till the size match
typedef SizedInteger<16>::Type int16_t;
typedef SizedUnsignedInteger<16>::Type uint16_t;
typedef SizedInteger<32>::Type int32_t;
typedef SizedUnsignedInteger<32>::Type uint32_t;
typedef SizedInteger<64>::Type int64_t;
typedef SizedUnsignedInteger<64>::Type uint64_t;
}
#endif

1
source/shared_lib/include/streflop/VERSION.txt Normal file
View File

@@ -0,0 +1 @@
0.2

196
source/shared_lib/include/streflop/X87DenormalSquasher.h Normal file
View File

@@ -0,0 +1,196 @@
/*
streflop: STandalone REproducible FLOating-Point
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
Heavily relies on GNU Libm, itself depending on netlib fplibm, GNU MP, and IBM MP lib.
Uses SoftFloat too.
Please read the history and copyright information in the documentation provided with the source code
*/
#ifndef X87DenormalSquasher_H
#define X87DenormalSquasher_H
/// This file should be included from within a streflop namespace
template<typename ValueType> inline void X87DenormalSquashFunction(ValueType& value) {
struct X {};
X Unknown_numeric_type;
// unknown types do not compile
value = Unknown_numeric_type;
}
// Use cmov instruction to avoid branching?
// Pro: branching may be costly, though the branching predictor may compensate when there are few denormals
// Con: cmov forces an unconditional writeback to the mem just after read, which may be worse than the branch
template<> inline void X87DenormalSquashFunction<float>(float& value) {
if ((reinterpret_cast<int*>(&value)[0] & 0x7F800000) == 0) value = 0.0f;
}
template<> inline void X87DenormalSquashFunction<double>(double& value) {
if ((reinterpret_cast<int*>(&value)[1] & 0x7FF00000) == 0) value = 0.0;
}
template<> inline void X87DenormalSquashFunction<long double>(long double& value) {
if ((reinterpret_cast<short*>(&value)[4] & 0x7FFF) == 0) value = 0.0L;
}
/// Wrapper class for the denormal squashing of X87
/// The goal is to be able to use this file for normal arithmetic operations
/// as if this was a native float or double, thanks to C++ operator overloading.
/// In particular, the GNU libm may be compiled this way, though no effort was made
/// to optimize out the denormal parts (and squashing denormals is non-IEEE compliant anyway)
template<typename ValueType> struct X87DenormalSquasher {
// Define value
ValueType value;
/// Uninitialized object
inline X87DenormalSquasher() {}
/// Copy and conversion from other precisions
inline X87DenormalSquasher(const X87DenormalSquasher<float>& f) : value(f.value) {X87DenormalSquashFunction<ValueType>(value);}
inline X87DenormalSquasher(const X87DenormalSquasher<double>& f) : value(f.value) {X87DenormalSquashFunction<ValueType>(value);}
inline X87DenormalSquasher(const X87DenormalSquasher<long double>& f) : value(f.value) {X87DenormalSquashFunction<ValueType>(value);}
inline X87DenormalSquasher& operator=(const X87DenormalSquasher<float>& f) {value = f.value; X87DenormalSquashFunction<ValueType>(value); return *this;} // self-ref OK
inline X87DenormalSquasher& operator=(const X87DenormalSquasher<double>& f) {value = f.value; X87DenormalSquashFunction<ValueType>(value); return *this;} // self-ref OK
inline X87DenormalSquasher& operator=(const X87DenormalSquasher<long double>& f) {value = f.value; X87DenormalSquashFunction<ValueType>(value); return *this;} // self-ref OK
/// Destructor
inline ~X87DenormalSquasher() {}
/// Now the real fun, arithmetic operator overloading
inline X87DenormalSquasher& operator+=(const X87DenormalSquasher& f) {value += f.value; X87DenormalSquashFunction<ValueType>(value); return *this;}
inline X87DenormalSquasher& operator-=(const X87DenormalSquasher& f) {value -= f.value; X87DenormalSquashFunction<ValueType>(value); return *this;}
inline X87DenormalSquasher& operator*=(const X87DenormalSquasher& f) {value *= f.value; X87DenormalSquashFunction<ValueType>(value); return *this;}
inline X87DenormalSquasher& operator/=(const X87DenormalSquasher& f) {value /= f.value; X87DenormalSquashFunction<ValueType>(value); return *this;}
inline bool operator==(const X87DenormalSquasher& f) const {return value == f.value;}
inline bool operator!=(const X87DenormalSquasher& f) const {return value != f.value;}
inline bool operator<(const X87DenormalSquasher& f) const {return value < f.value;}
inline bool operator<=(const X87DenormalSquasher& f) const {return value <= f.value;}
inline bool operator>(const X87DenormalSquasher& f) const {return value > f.value;}
inline bool operator>=(const X87DenormalSquasher& f) const {return value >= f.value;}
#define STREFLOP_X87DENORMAL_NATIVE_FLOAT_OPS(native_type) \
inline X87DenormalSquasher(const native_type f) {value = f; X87DenormalSquashFunction<ValueType>(value);} \
inline X87DenormalSquasher& operator=(const native_type f) {value = f; X87DenormalSquashFunction<ValueType>(value); return *this;} \
inline operator native_type() const {return static_cast<native_type>(value);} \
inline X87DenormalSquasher& operator+=(const native_type f) {return operator+=(X87DenormalSquasher(f));} \
inline X87DenormalSquasher& operator-=(const native_type f) {return operator-=(X87DenormalSquasher(f));} \
inline X87DenormalSquasher& operator*=(const native_type f) {return operator*=(X87DenormalSquasher(f));} \
inline X87DenormalSquasher& operator/=(const native_type f) {return operator/=(X87DenormalSquasher(f));} \
inline bool operator==(const native_type f) const {return operator==(X87DenormalSquasher(f));} \
inline bool operator!=(const native_type f) const {return operator!=(X87DenormalSquasher(f));} \
inline bool operator<(const native_type f) const {return operator<(X87DenormalSquasher(f));} \
inline bool operator<=(const native_type f) const {return operator<=(X87DenormalSquasher(f));} \
inline bool operator>(const native_type f) const {return operator>(X87DenormalSquasher(f));} \
inline bool operator>=(const native_type f) const {return operator>=(X87DenormalSquasher(f));}
#define STREFLOP_X87DENORMAL_NATIVE_INT_OPS(native_type) \
inline X87DenormalSquasher(const native_type x) {value = x;} \
inline X87DenormalSquasher& operator=(const native_type x) {value = x; return *this;} \
inline operator native_type() const {return static_cast<native_type>(value);} \
inline X87DenormalSquasher& operator+=(const native_type f) {value += f; X87DenormalSquashFunction<ValueType>(value); return *this;} \
inline X87DenormalSquasher& operator-=(const native_type f) {value -= f; X87DenormalSquashFunction<ValueType>(value); return *this;} \
inline X87DenormalSquasher& operator*=(const native_type f) {value *= f; X87DenormalSquashFunction<ValueType>(value); return *this;} \
inline X87DenormalSquasher& operator/=(const native_type f) {value /= f; X87DenormalSquashFunction<ValueType>(value); return *this;} \
inline bool operator==(const native_type f) const {return value == f;} \
inline bool operator!=(const native_type f) const {return value != f;} \
inline bool operator<(const native_type f) const {return value < f;} \
inline bool operator<=(const native_type f) const {return value <= f;} \
inline bool operator>(const native_type f) const {return value > f;} \
inline bool operator>=(const native_type f) const {return value >= f;}
STREFLOP_X87DENORMAL_NATIVE_FLOAT_OPS(float)
STREFLOP_X87DENORMAL_NATIVE_FLOAT_OPS(double)
STREFLOP_X87DENORMAL_NATIVE_FLOAT_OPS(long double)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(char)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(unsigned char)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(short)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(unsigned short)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(int)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(unsigned int)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(long)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(unsigned long)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(long long)
STREFLOP_X87DENORMAL_NATIVE_INT_OPS(unsigned long long)
};
/// binary operators
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator+(const X87DenormalSquasher<ValueType>& f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1.value + f2.value);}
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator-(const X87DenormalSquasher<ValueType>& f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1.value - f2.value);}
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator*(const X87DenormalSquasher<ValueType>& f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1.value * f2.value);}
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator/(const X87DenormalSquasher<ValueType>& f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1.value / f2.value);}
#define STREFLOP_X87DENORMAL_FLOAT_BINARY(native_type) \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator+(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return f1 + X87DenormalSquasher<ValueType>(f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator-(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return f1 - X87DenormalSquasher<ValueType>(f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator*(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return f1 * X87DenormalSquasher<ValueType>(f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator/(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return f1 / X87DenormalSquasher<ValueType>(f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator+(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1) + f2;} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator-(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1) - f2;} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator*(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1) * f2;} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator/(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1) / f2;} \
template<typename ValueType> inline bool operator==(const native_type value, const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(value) == f;} \
template<typename ValueType> inline bool operator!=(const native_type value, const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(value) != f;} \
template<typename ValueType> inline bool operator<(const native_type value, const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(value) < f;} \
template<typename ValueType> inline bool operator<=(const native_type value, const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(value) <= f;} \
template<typename ValueType> inline bool operator>(const native_type value, const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(value) > f;} \
template<typename ValueType> inline bool operator>=(const native_type value, const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(value) >= f;}
#define STREFLOP_X87DENORMAL_INT_BINARY(native_type) \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator+(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return X87DenormalSquasher<ValueType>(f1.value + f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator-(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return X87DenormalSquasher<ValueType>(f1.value - f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator*(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return X87DenormalSquasher<ValueType>(f1.value * f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator/(const X87DenormalSquasher<ValueType>& f1, const native_type f2) {return X87DenormalSquasher<ValueType>(f1.value / f2);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator+(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1 + f2.value);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator-(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1 - f2.value);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator*(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1 * f2.value);} \
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator/(const native_type f1, const X87DenormalSquasher<ValueType>& f2) {return X87DenormalSquasher<ValueType>(f1 / f2.value);} \
template<typename ValueType> inline bool operator==(const native_type value, const X87DenormalSquasher<ValueType>& f) {return value == f.value;} \
template<typename ValueType> inline bool operator!=(const native_type value, const X87DenormalSquasher<ValueType>& f) {return value != f.value;} \
template<typename ValueType> inline bool operator<(const native_type value, const X87DenormalSquasher<ValueType>& f) {return value < f.value;} \
template<typename ValueType> inline bool operator<=(const native_type value, const X87DenormalSquasher<ValueType>& f) {return value <= f.value;} \
template<typename ValueType> inline bool operator>(const native_type value, const X87DenormalSquasher<ValueType>& f) {return value > f.value;} \
template<typename ValueType> inline bool operator>=(const native_type value, const X87DenormalSquasher<ValueType>& f) {return value >= f.value;}
STREFLOP_X87DENORMAL_FLOAT_BINARY(float)
STREFLOP_X87DENORMAL_FLOAT_BINARY(double)
STREFLOP_X87DENORMAL_FLOAT_BINARY(long double)
STREFLOP_X87DENORMAL_INT_BINARY(char)
STREFLOP_X87DENORMAL_INT_BINARY(unsigned char)
STREFLOP_X87DENORMAL_INT_BINARY(short)
STREFLOP_X87DENORMAL_INT_BINARY(unsigned short)
STREFLOP_X87DENORMAL_INT_BINARY(int)
STREFLOP_X87DENORMAL_INT_BINARY(unsigned int)
STREFLOP_X87DENORMAL_INT_BINARY(long)
STREFLOP_X87DENORMAL_INT_BINARY(unsigned long)
STREFLOP_X87DENORMAL_INT_BINARY(long long)
STREFLOP_X87DENORMAL_INT_BINARY(unsigned long long)
/// Unary operators
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator-(const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(-f.value);}
template<typename ValueType> inline X87DenormalSquasher<ValueType> operator+(const X87DenormalSquasher<ValueType>& f) {return X87DenormalSquasher<ValueType>(+f.value);}
/*
/// Stream operators for output
/// Warning! The SYSTEM operators are used to the normal floats and therefore
/// this function may not be reliable (they may depend on system rounding mode, etc.)
/// Still useful for debug output and to report results with a few significant digits
template<typename ValueType> inline std::ostream& operator<<(std::ostream& out, const X87DenormalSquasher<ValueType>& f) {
out << f.value;
return out;
}
template<typename ValueType> inline std::istream& operator>>(std::istream& in, X87DenormalSquasher<ValueType>& f) {
in >> f.value;
return in;
}
*/
#endif

19
source/shared_lib/include/streflop/libm/README.txt Normal file
View File

@@ -0,0 +1,19 @@
This directory contains (in addition to this file):
- flt-32/ dbl-64/ ldbl-96/ headers/: These subdirectories are imported from the GNU libm automatically. Do not modify them manually. You may safely erase them and run the import.pl script again. The "clean-import" make target does just that.
- import.pl: A perl script that takes as argument the location of a valid glibc-2.4 source directory. It will create the aforementioned subdirectories and convert the original libm files. The conversion mainly concerns removing external dependencies, and making the c files and constructs compilable with the streflop wrapper types. See the import.pl script itself for the full list of modifications.
- Makefile: Instructions for compiling the subdirectories with the make command. You you re-use this as the basis for other build systems (like CMake, scons, etc.).
- streflop_libm_bridge.h: Contains the definitions and other macros that are necessary to compile libm, in a streflop framework. These declarations were either done by external files, or by other parts of the glibc that were not imported.
- e_expf.c: The float exp in e_expf.c uses doubles internally since revision 1.2 in the libm-ieee754 CVS attic! This is the slower, but purely float version, that is rolled back from the CVS attic.
- w_expf.c: A wrapper to expf that replaces the libm wrapper by a wraper to the float only version.
- (after compilation): flt-target dbl-target ldbl-target temporary files for the make process
The original GNU libm is released under the GNU LGPL license, and so are these modifications. See the LGPL.txt in the parent streflop main directory. See also the comments at the beginning of each file for particular information, especially the Sun Microsystems disclaimer.
If you import a version of glibc other than 2.4, please check the information there (in particular their LICENSES file) for potential inclusion of new files in the libm subdirectories that are used here.

5
source/shared_lib/include/streflop/libm/flt-32/Makefile Normal file
View File

@@ -0,0 +1,5 @@
# Makefile automatically generated by import.pl
include ../../Makefile.common
CPPFLAGS += -I../headers -DLIBM_COMPILING_FLT32=1
all: e_acosf.o e_acoshf.o e_asinf.o e_atan2f.o e_atanhf.o e_coshf.o e_exp2f.o e_expf.o e_fmodf.o e_gammaf_r.o e_hypotf.o e_j0f.o e_j1f.o e_jnf.o e_lgammaf_r.o e_log10f.o e_log2f.o e_logf.o e_powf.o e_rem_pio2f.o e_remainderf.o e_sinhf.o e_sqrtf.o k_cosf.o k_rem_pio2f.o k_sinf.o k_tanf.o s_asinhf.o s_atanf.o s_cbrtf.o s_ceilf.o s_copysignf.o s_cosf.o s_erff.o s_expm1f.o s_fabsf.o s_finitef.o s_floorf.o s_fpclassifyf.o s_frexpf.o s_ilogbf.o s_isinff.o s_isnanf.o s_ldexpf.o s_llrintf.o s_llroundf.o s_log1pf.o s_logbf.o s_lrintf.o s_lroundf.o s_modff.o s_nearbyintf.o s_nextafterf.o s_remquof.o s_rintf.o s_roundf.o s_scalblnf.o s_scalbnf.o s_signbitf.o s_sincosf.o s_sinf.o s_tanf.o s_tanhf.o s_truncf.o w_expf.o
echo 'flt-32 done!'

355
source/shared_lib/include/streflop/libm/flt-32/t_exp2f.h Normal file
View File

@@ -0,0 +1,355 @@
/* See the import.pl script for potential modifications */
/* Accurate tables for exp2f().
Copyright (C) 1998 Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Geoffrey Keating <geoffk@ozemail.com.au>
The GNU C Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1f of the License, or (at your option) any later version.
The GNU C Library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with the GNU C Library; if not, write to the Free
Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
02111-1307 USA. */
/* This table has the property that, for all integers -128 <= i <= 127,
exp(i/256.0f + __exp2f_deltatable[i-128]) == __exp2f_atable[i+128] + r
for some -2^-35 < r < 2^-35 (abs(r) < 2^-36 if i <= 0); and that
__exp2f_deltatable[i+128] == t * 2^-30
for integer t so that abs(t) <= 43447 * 2^0. */
#define W30 (9.31322575e-10f)
static const Simple __exp2f_deltatable[256] = {
-810*W30, 283*W30, -1514*W30, 1304*W30,
-1148*W30, -98*W30, -744*W30, -156*W30,
-419*W30, -155*W30, 474*W30, 167*W30,
-1984*W30, -826*W30, 692*W30, 781*W30,
-578*W30, -411*W30, -129*W30, -1500*W30,
654*W30, -141*W30, -816*W30, -53*W30,
148*W30, 493*W30, -2214*W30, 760*W30,
260*W30, 750*W30, -1300*W30, 1424*W30,
-1445*W30, -339*W30, -680*W30, -349*W30,
-922*W30, 531*W30, 193*W30, -2892*W30,
290*W30, -2145*W30, -276*W30, 485*W30,
-695*W30, 215*W30, -7093*W30, 412*W30,
-4596*W30, 367*W30, 592*W30, -615*W30,
-97*W30, -1066*W30, 972*W30, -226*W30,
-625*W30, -374*W30, -5647*W30, -180*W30,
20349*W30, -447*W30, 111*W30, -4164*W30,
-87*W30, -21*W30, -251*W30, 66*W30,
-517*W30, 2093*W30, -263*W30, 182*W30,
-601*W30, 475*W30, -483*W30, -1251*W30,
-373*W30, 1471*W30, -92*W30, -215*W30,
-97*W30, -190*W30, 0*W30, -290*W30,
-2647*W30, 1940*W30, -582*W30, 28*W30,
833*W30, 1493*W30, 34*W30, 321*W30,
3327*W30, -35*W30, 177*W30, -135*W30,
-796*W30, -428*W30, 129*W30, 9332*W30,
-12*W30, -69*W30, -1743*W30, 6508*W30,
-60*W30, 359*W30, 43447*W30, 15*W30,
-23*W30, -305*W30, -375*W30, -652*W30,
667*W30, 269*W30, -1575*W30, 185*W30,
-329*W30, 200*W30, 6002*W30, 163*W30,
-647*W30, 19*W30, -603*W30, -755*W30,
742*W30, -438*W30, 3587*W30, 2560*W30,
0*W30, -520*W30, -241*W30, -299*W30,
-1270*W30, -991*W30, -1138*W30, 255*W30,
-1192*W30, 1722*W30, 1023*W30, 3700*W30,
-1388*W30, -1551*W30, -2549*W30, 27*W30,
282*W30, 673*W30, 113*W30, 1561*W30,
72*W30, 873*W30, 87*W30, -395*W30,
-433*W30, 629*W30, 3440*W30, -284*W30,
-592*W30, -103*W30, -46*W30, -3844*W30,
1712*W30, 303*W30, 1555*W30, -631*W30,
-1400*W30, -961*W30, -854*W30, -276*W30,
407*W30, 833*W30, -345*W30, -1501*W30,
121*W30, -1581*W30, 400*W30, 150*W30,
1224*W30, -139*W30, -563*W30, 879*W30,
933*W30, 2939*W30, 788*W30, 211*W30,
530*W30, -192*W30, 706*W30, -13347*W30,
1065*W30, 3*W30, 111*W30, -208*W30,
-360*W30, -532*W30, -291*W30, 483*W30,
987*W30, -33*W30, -1373*W30, -166*W30,
-1174*W30, -3955*W30, 1601*W30, -280*W30,
1405*W30, 600*W30, -1659*W30, -23*W30,
390*W30, 449*W30, 570*W30, -13143*W30,
-9*W30, -1646*W30, 1201*W30, 294*W30,
2181*W30, -1173*W30, 1388*W30, -4504*W30,
190*W30, -2304*W30, 211*W30, 239*W30,
48*W30, -817*W30, 1018*W30, 1828*W30,
-663*W30, 1408*W30, 408*W30, -36*W30,
1295*W30, -230*W30, 1341*W30, 9*W30,
40*W30, 705*W30, 186*W30, 376*W30,
557*W30, 5866*W30, 363*W30, -1558*W30,
718*W30, 669*W30, 1369*W30, -2972*W30,
-468*W30, -121*W30, -219*W30, 667*W30,
29954*W30, 366*W30, 48*W30, -203*W30
};
namespace streflop_libm {
static const Simple __exp2f_atable[256] /* __attribute__((mode(SF))) */ = {
0.707106411447f, /* 0x0.b504ecfff */
0.709024071690f, /* 0x0.b58299fff */
0.710945606239f, /* 0x0.b60088000 */
0.712874472142f, /* 0x0.b67ef1000 */
0.714806139464f, /* 0x0.b6fd88fff */
0.716744661340f, /* 0x0.b77c94000 */
0.718687653549f, /* 0x0.b7fbea000 */
0.720636486992f, /* 0x0.b87ba1fff */
0.722590208040f, /* 0x0.b8fbabfff */
0.724549472323f, /* 0x0.b97c12fff */
0.726514220228f, /* 0x0.b9fcd5fff */
0.728483855735f, /* 0x0.ba7deb000 */
0.730457961549f, /* 0x0.baff4afff */
0.732438981522f, /* 0x0.bb811efff */
0.734425544748f, /* 0x0.bc0350000 */
0.736416816713f, /* 0x0.bc85d0000 */
0.738412797450f, /* 0x0.bd089efff */
0.740414917465f, /* 0x0.bd8bd4fff */
0.742422521111f, /* 0x0.be0f66fff */
0.744434773914f, /* 0x0.be9346fff */
0.746454179287f, /* 0x0.bf179f000 */
0.748477637755f, /* 0x0.bf9c3afff */
0.750506639473f, /* 0x0.c02133fff */
0.752541840064f, /* 0x0.c0a694fff */
0.754582285889f, /* 0x0.c12c4e000 */
0.756628334525f, /* 0x0.c1b265000 */
0.758678436269f, /* 0x0.c238bffff */
0.760736882681f, /* 0x0.c2bfa6fff */
0.762799203401f, /* 0x0.c346cf000 */
0.764867603790f, /* 0x0.c3ce5d000 */
0.766940355298f, /* 0x0.c45633fff */
0.769021093841f, /* 0x0.c4de90fff */
0.771104693409f, /* 0x0.c5671dfff */
0.773195922364f, /* 0x0.c5f02afff */
0.775292098512f, /* 0x0.c6798afff */
0.777394294745f, /* 0x0.c70350000 */
0.779501736166f, /* 0x0.c78d6d000 */
0.781615912910f, /* 0x0.c817fafff */
0.783734917628f, /* 0x0.c8a2d9fff */
0.785858273516f, /* 0x0.c92e02000 */
0.787990570071f, /* 0x0.c9b9c0000 */
0.790125787245f, /* 0x0.ca45aefff */
0.792268991467f, /* 0x0.cad223fff */
0.794417440881f, /* 0x0.cb5ef0fff */
0.796570718287f, /* 0x0.cbec0efff */
0.798730909811f, /* 0x0.cc79a0fff */
0.800892710672f, /* 0x0.cd074dfff */
0.803068041795f, /* 0x0.cd95ddfff */
0.805242776881f, /* 0x0.ce2464000 */
0.807428598393f, /* 0x0.ceb3a3fff */
0.809617877002f, /* 0x0.cf431dfff */
0.811812341211f, /* 0x0.cfd2eefff */
0.814013659956f, /* 0x0.d06333000 */
0.816220164311f, /* 0x0.d0f3ce000 */
0.818434238424f, /* 0x0.d184e7fff */
0.820652604094f, /* 0x0.d21649fff */
0.822877407074f, /* 0x0.d2a818000 */
0.825108587751f, /* 0x0.d33a51000 */
0.827342867839f, /* 0x0.d3ccbdfff */
0.829588949684f, /* 0x0.d45ff1000 */
0.831849217401f, /* 0x0.d4f411fff */
0.834093391880f, /* 0x0.d58724fff */
0.836355149750f, /* 0x0.d61b5f000 */
0.838620424257f, /* 0x0.d6afd3fff */
0.840896368027f, /* 0x0.d744fc000 */
0.843176305293f, /* 0x0.d7da66fff */
0.845462262643f, /* 0x0.d87037000 */
0.847754716864f, /* 0x0.d90673fff */
0.850052893157f, /* 0x0.d99d10fff */
0.852359056469f, /* 0x0.da3433fff */
0.854668736446f, /* 0x0.dacb91fff */
0.856986224651f, /* 0x0.db6373000 */
0.859309315673f, /* 0x0.dbfbb1fff */
0.861639738080f, /* 0x0.dc946bfff */
0.863975346095f, /* 0x0.dd2d7d000 */
0.866317391394f, /* 0x0.ddc6f9fff */
0.868666708472f, /* 0x0.de60f1000 */
0.871022939695f, /* 0x0.defb5c000 */
0.873383641229f, /* 0x0.df9611fff */
0.875751554968f, /* 0x0.e03141000 */
0.878126025200f, /* 0x0.e0ccde000 */
0.880506813521f, /* 0x0.e168e4fff */
0.882894217966f, /* 0x0.e2055afff */
0.885287821299f, /* 0x0.e2a239000 */
0.887686729423f, /* 0x0.e33f6ffff */
0.890096127973f, /* 0x0.e3dd56fff */
0.892507970338f, /* 0x0.e47b67000 */
0.894928157336f, /* 0x0.e51a03000 */
0.897355020043f, /* 0x0.e5b90efff */
0.899788379682f, /* 0x0.e65888000 */
0.902227103705f, /* 0x0.e6f85afff */
0.904673457151f, /* 0x0.e798ae000 */
0.907128036008f, /* 0x0.e8398afff */
0.909585535528f, /* 0x0.e8da99000 */
0.912051796915f, /* 0x0.e97c3a000 */
0.914524436003f, /* 0x0.ea1e46000 */
0.917003571999f, /* 0x0.eac0bf000 */
0.919490039339f, /* 0x0.eb63b2fff */
0.921983361257f, /* 0x0.ec071a000 */
0.924488604054f, /* 0x0.ecab48fff */
0.926989555360f, /* 0x0.ed4f30000 */
0.929502844812f, /* 0x0.edf3e6000 */
0.932021975503f, /* 0x0.ee98fdfff */
0.934553921208f, /* 0x0.ef3eecfff */
0.937083780759f, /* 0x0.efe4b8fff */
0.939624726786f, /* 0x0.f08b3f000 */
0.942198514924f, /* 0x0.f133ebfff */
0.944726586343f, /* 0x0.f1d99a000 */
0.947287976728f, /* 0x0.f28176fff */
0.949856162070f, /* 0x0.f329c5fff */
0.952431440345f, /* 0x0.f3d28bfff */
0.955013573175f, /* 0x0.f47bc5000 */
0.957603693021f, /* 0x0.f52584000 */
0.960199773321f, /* 0x0.f5cfa7000 */
0.962801992906f, /* 0x0.f67a31000 */
0.965413510788f, /* 0x0.f72556fff */
0.968030691152f, /* 0x0.f7d0dc000 */
0.970655620084f, /* 0x0.f87ce2fff */
0.973290979849f, /* 0x0.f92998fff */
0.975926160805f, /* 0x0.f9d64bfff */
0.978571653370f, /* 0x0.fa83ac000 */
0.981225252139f, /* 0x0.fb3193fff */
0.983885228626f, /* 0x0.fbdfe6fff */
0.986552715296f, /* 0x0.fc8eb7fff */
0.989228487027f, /* 0x0.fd3e14000 */
0.991909801964f, /* 0x0.fdedcd000 */
0.994601726545f, /* 0x0.fe9e38000 */
0.997297704209f, /* 0x0.ff4ee6fff */
1.000000000000f, /* 0x1.000000000f */
1.002710938457f, /* 0x1.00b1aa000 */
1.005429744692f, /* 0x1.0163d7ffe */
1.008155703526f, /* 0x1.02167dffe */
1.010888457284f, /* 0x1.02c995fff */
1.013629436498f, /* 0x1.037d38000 */
1.016377568250f, /* 0x1.043152000f */
1.019134163841f, /* 0x1.04e5f9ffe */
1.021896362316f, /* 0x1.059b00000 */
1.024668931945f, /* 0x1.0650b3ffe */
1.027446627635f, /* 0x1.0706be001 */
1.030234098408f, /* 0x1.07bd6bffe */
1.033023953416f, /* 0x1.087441ffe */
1.035824656494f, /* 0x1.092bce000 */
1.038632392900f, /* 0x1.09e3d0001 */
1.041450142840f, /* 0x1.0a9c79ffe */
1.044273972530f, /* 0x1.0b558a001 */
1.047105550795f, /* 0x1.0c0f1c001 */
1.049944162390f, /* 0x1.0cc924001 */
1.052791833895f, /* 0x1.0d83c4001 */
1.055645227426f, /* 0x1.0e3ec3fff */
1.058507919326f, /* 0x1.0efa60001 */
1.061377286898f, /* 0x1.0fb66bfff */
1.064254641510f, /* 0x1.1072fdffe */
1.067140102389f, /* 0x1.113018000f */
1.070034146304f, /* 0x1.11edc1fff */
1.072937250162f, /* 0x1.12ac04001 */
1.075843691823f, /* 0x1.136a7dfff */
1.078760385496f, /* 0x1.1429a3ffe */
1.081685543070f, /* 0x1.14e958000f */
1.084618330005f, /* 0x1.15a98c000 */
1.087556362176f, /* 0x1.166a18001 */
1.090508937863f, /* 0x1.172b98001 */
1.093464612954f, /* 0x1.17ed4bfff */
1.096430182434f, /* 0x1.18afa5ffe */
1.099401354802f, /* 0x1.19725e000f */
1.102381587017f, /* 0x1.1a35adfff */
1.105370759965f, /* 0x1.1af994000 */
1.108367800686f, /* 0x1.1bbdfdffe */
1.111373305331f, /* 0x1.1c82f6000 */
1.114387035385f, /* 0x1.1d4878001 */
1.117408752440f, /* 0x1.1e0e7ffff */
1.120437502874f, /* 0x1.1ed4fe000 */
1.123474478729f, /* 0x1.1f9c06000 */
1.126521706601f, /* 0x1.2063ba001 */
1.129574775716f, /* 0x1.212bd0001 */
1.132638812065f, /* 0x1.21f49e000 */
1.135709524130f, /* 0x1.22bddbffe */
1.138789534565f, /* 0x1.2387b5fff */
1.141876101508f, /* 0x1.2451fe000 */
1.144971728301f, /* 0x1.251cddffe */
1.148077130296f, /* 0x1.25e861ffe */
1.151189923305f, /* 0x1.26b462001 */
1.154312610610f, /* 0x1.278107ffe */
1.157440662410f, /* 0x1.284e08001f */
1.160578370109f, /* 0x1.291baa001 */
1.163725256932f, /* 0x1.29e9e6000 */
1.166879892324f, /* 0x1.2ab8a3ffe */
1.170044302935f, /* 0x1.2b8805fff */
1.173205971694f, /* 0x1.2c5739ffe */
1.176397800428f, /* 0x1.2d2867ffe */
1.179586529747f, /* 0x1.2df962001 */
1.182784795737f, /* 0x1.2ecafbffe */
1.185991406414f, /* 0x1.2f9d21ffe */
1.189206838636f, /* 0x1.306fdc001 */
1.192430973067f, /* 0x1.314328000f */
1.195664167430f, /* 0x1.32170c001 */
1.198906540890f, /* 0x1.32eb8a001 */
1.202157497408f, /* 0x1.33c098000 */
1.205416083326f, /* 0x1.349625fff */
1.208683252332f, /* 0x1.356c43fff */
1.211961269402f, /* 0x1.364318001f */
1.215246438983f, /* 0x1.371a64000 */
1.218539118740f, /* 0x1.37f22dffe */
1.221847295770f, /* 0x1.38cafc000 */
1.225158572187f, /* 0x1.39a3fdfff */
1.228481650325f, /* 0x1.3a7dc5ffe */
1.231811761846f, /* 0x1.3b5803fff */
1.235149741144f, /* 0x1.3c32c5ffe */
1.238499879811f, /* 0x1.3d0e53ffe */
1.241858124726f, /* 0x1.3dea69fff */
1.245225191102f, /* 0x1.3ec713fff */
1.248601436624f, /* 0x1.3fa458000 */
1.251975655584f, /* 0x1.40817a001 */
1.255380749731f, /* 0x1.4160a2001 */
1.258783102010f, /* 0x1.423f9bffe */
1.262198328973f, /* 0x1.431f6e000 */
1.265619754780f, /* 0x1.43ffa7fff */
1.269052743928f, /* 0x1.44e0a4001 */
1.272490739830f, /* 0x1.45c1f4000 */
1.275942921659f, /* 0x1.46a432001 */
1.279397487615f, /* 0x1.478697ffe */
1.282870173427f, /* 0x1.486a2dffe */
1.286346316319f, /* 0x1.494dfdffe */
1.289836049094f, /* 0x1.4a32b2001 */
1.293333172770f, /* 0x1.4b17e1ffe */
1.296839594835f, /* 0x1.4bfdadfff */
1.300354957560f, /* 0x1.4ce40fffe */
1.303882122055f, /* 0x1.4dcb38001 */
1.307417988757f, /* 0x1.4eb2f1ffe */
1.310960650439f, /* 0x1.4f9b1dfff */
1.314516782746f, /* 0x1.50842bfff */
1.318079948424f, /* 0x1.516daffff */
1.321653246888f, /* 0x1.5257de000 */
1.325237751030f, /* 0x1.5342c8001 */
1.328829526907f, /* 0x1.542e2c000 */
1.332433700535f, /* 0x1.551a5fffe */
1.336045145966f, /* 0x1.56070dffe */
1.339667558645f, /* 0x1.56f473ffe */
1.343300342533f, /* 0x1.57e287ffe */
1.346941947961f, /* 0x1.58d130001 */
1.350594043714f, /* 0x1.59c087ffe */
1.354256033883f, /* 0x1.5ab085fff */
1.357932448365f, /* 0x1.5ba175ffe */
1.361609339707f, /* 0x1.5c926dfff */
1.365299344044f, /* 0x1.5d8441ffe */
1.369003057507f, /* 0x1.5e76fc001 */
1.372714757920f, /* 0x1.5f6a3c000 */
1.376437187179f, /* 0x1.605e2fffe */
1.380165219333f, /* 0x1.615282001f */
1.383909463864f, /* 0x1.6247e3ffe */
1.387661933907f, /* 0x1.633dd0000 */
1.391424179060f, /* 0x1.64345fffe */
1.395197510706f, /* 0x1.652ba9fff */
1.399006724329f, /* 0x1.66254dffe */
1.402773022651f, /* 0x1.671c22000 */
1.406576037403f, /* 0x1.68155dfff */
1.410389423392f, /* 0x1.690f48001 */
};
}

7
source/shared_lib/include/streflop/libm/headers/SMath.h Normal file
View File

@@ -0,0 +1,7 @@
/* Copied from math.h in same directory */
/* This file is stub automatically generated by import.pl */
// Include bridge, just in case math_private is not included
#include "../streflop_libm_bridge.h"

77
source/shared_lib/include/streflop/libm/headers/endian.h Normal file
View File

@@ -0,0 +1,77 @@
/* See the import.pl script for potential modifications */
/* Copyright (C) 1992, 1996, 1997, 2000 Free Software Foundation, Inc.
This file is part of the GNU C Library.
The GNU C Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1 of the License, or (at your option) any later version.
The GNU C Library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with the GNU C Library; if not, write to the Free
Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
02111-1307 USA. */
#ifndef _ENDIAN_H
#define _ENDIAN_H 1
#include "features.h"
/* Definitions for byte order, according to significance of bytes,
from low addresses to high addresses. The value is what you get by
putting '4' in the most significant byte, '3' in the second most
significant byte, '2' in the second least significant byte, and '1'
in the least significant byte, and then writing down one digit for
each byte, starting with the byte at the lowest address at the left,
and proceeding to the byte with the highest address at the right. */
#define __LITTLE_ENDIAN 1234
#define __BIG_ENDIAN 4321
#define __PDP_ENDIAN 3412
/* This file defines `__BYTE_ORDER' for the particular machine. */
#include "../streflop_libm_bridge.h"
/* Some machines may need to use a different endianness for floating point
values. */
#ifndef __FLOAT_WORD_ORDER
# define __FLOAT_WORD_ORDER __BYTE_ORDER
#endif
#ifdef __USE_BSD
# define LITTLE_ENDIAN __LITTLE_ENDIAN
# define BIG_ENDIAN __BIG_ENDIAN
# define PDP_ENDIAN __PDP_ENDIAN
# define BYTE_ORDER __BYTE_ORDER
#endif
#if __BYTE_ORDER == __LITTLE_ENDIAN
# define __LONG_LONG_PAIR(HI, LO) LO, HI
#elif __BYTE_ORDER == __BIG_ENDIAN
# define __LONG_LONG_PAIR(HI, LO) HI, LO
#endif
#endif /* endian.h */
//#include <string/endian.h>
#if 1
//#if defined _LIBC && !defined _ISOMAC
# if __FLOAT_WORD_ORDER == __BIG_ENDIAN
# define BIG_ENDI 1
# undef LITTLE_ENDI
# define HIGH_HALF 0
# define LOW_HALF 1
# else
# if __FLOAT_WORD_ORDER == __LITTLE_ENDIAN
# undef BIG_ENDI
# define LITTLE_ENDI 1
# define HIGH_HALF 1
# define LOW_HALF 0
# endif
# endif
#endif

350
source/shared_lib/include/streflop/libm/headers/features.h Normal file
View File

@@ -0,0 +1,350 @@
/* See the import.pl script for potential modifications */
/* Copyright (C) 1991,1992,1993,1995,1996,1997,1998,1999,2000,2001,2002,2003,2004,2005,2006
Free Software Foundation, Inc.
This file is part of the GNU C Library.
The GNU C Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1 of the License, or (at your option) any later version.
The GNU C Library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with the GNU C Library; if not, write to the Free
Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
02111-1307 USA. */
#ifndef _FEATURES_H
#define _FEATURES_H 1
/* These are defined by the user (or the compiler)
to specify the desired environment:
__STRICT_ANSI__ ISO Standard C.
_ISOC99_SOURCE Extensions to ISO C89 from ISO C99.
_POSIX_SOURCE IEEE Std 1003.1.
_POSIX_C_SOURCE If ==1, like _POSIX_SOURCE; if >=2 add IEEE Std 1003.2;
if >=199309L, add IEEE Std 1003.1b-1993;
if >=199506L, add IEEE Std 1003.1c-1995;
if >=200112L, all of IEEE 1003.1-2004
_XOPEN_SOURCE Includes POSIX and XPG things. Set to 500 if
Single Unix conformance is wanted, to 600 for the
upcoming sixth revision.
_XOPEN_SOURCE_EXTENDED XPG things and X/Open Unix extensions.
_LARGEFILE_SOURCE Some more functions for correct standard I/O.
_LARGEFILE64_SOURCE Additional functionality from LFS for large files.
_FILE_OFFSET_BITS=N Select default filesystem interface.
_BSD_SOURCE ISO C, POSIX, and 4.3BSD things.
_SVID_SOURCE ISO C, POSIX, and SVID things.
_ATFILE_SOURCE Additional *at interfaces.
_GNU_SOURCE All of the above, plus GNU extensions.
_REENTRANT Select additionally reentrant object.
_THREAD_SAFE Same as _REENTRANT, often used by other systems.
_FORTIFY_SOURCE If set to numeric value > 0 additional security
measures are defined, according to level.
The `-ansi' switch to the GNU C compiler defines __STRICT_ANSI__.
If none of these are defined, the default is to have _SVID_SOURCE,
_BSD_SOURCE, and _POSIX_SOURCE set to one and _POSIX_C_SOURCE set to
199506L. If more than one of these are defined, they accumulate.
For example __STRICT_ANSI__, _POSIX_SOURCE and _POSIX_C_SOURCE
together give you ISO C, 1003.1, and 1003.2, but nothing else.
These are defined by this file and are used by the
header files to decide what to declare or define:
__USE_ISOC99 Define ISO C99 things.
__USE_POSIX Define IEEE Std 1003.1 things.
__USE_POSIX2 Define IEEE Std 1003.2 things.
__USE_POSIX199309 Define IEEE Std 1003.1, and .1b things.
__USE_POSIX199506 Define IEEE Std 1003.1, .1b, .1c and .1i things.
__USE_XOPEN Define XPG things.
__USE_XOPEN_EXTENDED Define X/Open Unix things.
__USE_UNIX98 Define Single Unix V2 things.
__USE_XOPEN2K Define XPG6 things.
__USE_LARGEFILE Define correct standard I/O things.
__USE_LARGEFILE64 Define LFS things with separate names.
__USE_FILE_OFFSET64 Define 64bit interface as default.
__USE_BSD Define 4.3BSD things.
__USE_SVID Define SVID things.
__USE_MISC Define things common to BSD and System V Unix.
__USE_ATFILE Define *at interfaces and AT_* constants for them.
__USE_GNU Define GNU extensions.
__USE_REENTRANT Define reentrant/thread-safe *_r functions.
__USE_FORTIFY_LEVEL Additional security measures used, according to level.
__FAVOR_BSD Favor 4.3BSD things in cases of conflict.
The macros `__GNU_LIBRARY__', `__GLIBC__', and `__GLIBC_MINOR__' are
defined by this file unconditionally. `__GNU_LIBRARY__' is provided
only for compatibility. All new code should use the other symbols
to test for features.
All macros listed above as possibly being defined by this file are
explicitly undefined if they are not explicitly defined.
Feature-test macros that are not defined by the user or compiler
but are implied by the other feature-test macros defined (or by the
lack of any definitions) are defined by the file. */
/* Undefine everything, so we get a clean slate. */
#undef __USE_ISOC99
#undef __USE_POSIX
#undef __USE_POSIX2
#undef __USE_POSIX199309
#undef __USE_POSIX199506
#undef __USE_XOPEN
#undef __USE_XOPEN_EXTENDED
#undef __USE_UNIX98
#undef __USE_XOPEN2K
#undef __USE_LARGEFILE
#undef __USE_LARGEFILE64
#undef __USE_FILE_OFFSET64
#undef __USE_BSD
#undef __USE_SVID
#undef __USE_MISC
#undef __USE_ATFILE
#undef __USE_GNU
#undef __USE_REENTRANT
#undef __USE_FORTIFY_LEVEL
#undef __FAVOR_BSD
#undef __KERNEL_STRICT_NAMES
/* Suppress kernel-name space pollution unless user expressedly asks
for it. */
#ifndef _LOOSE_KERNEL_NAMES
# define __KERNEL_STRICT_NAMES
#endif
/* Always use ISO C things. */
#define __USE_ANSI 1
/* Convenience macros to test the versions of glibc and gcc.
Use them like this:
#if __GNUC_PREREQ (2,8)
... code requiring gcc 2.8 or later ...
#endif
Note - they won't work for gcc1 or glibc1, since the _MINOR macros
were not defined then. */
#if defined __GNUC__ && defined __GNUC_MINOR__
# define __GNUC_PREREQ(maj, min) \
((__GNUC__ << 16) + __GNUC_MINOR__ >= ((maj) << 16) + (min))
#else
# define __GNUC_PREREQ(maj, min) 0
#endif
/* If _BSD_SOURCE was defined by the user, favor BSD over POSIX. */
#if defined _BSD_SOURCE && \
!(defined _POSIX_SOURCE || defined _POSIX_C_SOURCE || \
defined _XOPEN_SOURCE || defined _XOPEN_SOURCE_EXTENDED || \
defined _GNU_SOURCE || defined _SVID_SOURCE)
# define __FAVOR_BSD 1
#endif
/* If _GNU_SOURCE was defined by the user, turn on all the other features. */
#ifdef _GNU_SOURCE
# undef _ISOC99_SOURCE
# define _ISOC99_SOURCE 1
# undef _POSIX_SOURCE
# define _POSIX_SOURCE 1
# undef _POSIX_C_SOURCE
# define _POSIX_C_SOURCE 199506L
# undef _XOPEN_SOURCE
# define _XOPEN_SOURCE 600
# undef _XOPEN_SOURCE_EXTENDED
# define _XOPEN_SOURCE_EXTENDED 1
# undef _LARGEFILE64_SOURCE
# define _LARGEFILE64_SOURCE 1
# undef _BSD_SOURCE
# define _BSD_SOURCE 1
# undef _SVID_SOURCE
# define _SVID_SOURCE 1
# undef _ATFILE_SOURCE
# define _ATFILE_SOURCE 1
#endif
/* If nothing (other than _GNU_SOURCE) is defined,
define _BSD_SOURCE and _SVID_SOURCE. */
#if (!defined __STRICT_ANSI__ && !defined _ISOC99_SOURCE && \
!defined _POSIX_SOURCE && !defined _POSIX_C_SOURCE && \
!defined _XOPEN_SOURCE && !defined _XOPEN_SOURCE_EXTENDED && \
!defined _BSD_SOURCE && !defined _SVID_SOURCE)
# define _BSD_SOURCE 1
# define _SVID_SOURCE 1
#endif
/* This is to enable the ISO C99 extension. Also recognize the old macro
which was used prior to the standard acceptance. This macro will
eventually go away and the features enabled by default once the ISO C99
standard is widely adopted. */
#if (defined _ISOC99_SOURCE || defined _ISOC9X_SOURCE \
|| (defined __STDC_VERSION__ && __STDC_VERSION__ >= 199901L))
# define __USE_ISOC99 1
#endif
/* If none of the ANSI/POSIX macros are defined, use POSIX.1 and POSIX.2
(and IEEE Std 1003.1b-1993 unless _XOPEN_SOURCE is defined). */
#if ((!defined __STRICT_ANSI__ || (_XOPEN_SOURCE - 0) >= 500) && \
!defined _POSIX_SOURCE && !defined _POSIX_C_SOURCE)
# define _POSIX_SOURCE 1
# if defined _XOPEN_SOURCE && (_XOPEN_SOURCE - 0) < 500
# define _POSIX_C_SOURCE 2
# elif defined _XOPEN_SOURCE && (_XOPEN_SOURCE - 0) < 600
# define _POSIX_C_SOURCE 199506L
# else
# define _POSIX_C_SOURCE 200112L
# endif
#endif
#if defined _POSIX_SOURCE || _POSIX_C_SOURCE >= 1 || defined _XOPEN_SOURCE
# define __USE_POSIX 1
#endif
#if defined _POSIX_C_SOURCE && _POSIX_C_SOURCE >= 2 || defined _XOPEN_SOURCE
# define __USE_POSIX2 1
#endif
#if (_POSIX_C_SOURCE - 0) >= 199309L
# define __USE_POSIX199309 1
#endif
#if (_POSIX_C_SOURCE - 0) >= 199506L
# define __USE_POSIX199506 1
#endif
#if (_POSIX_C_SOURCE - 0) >= 200112L
# define __USE_XOPEN2K 1
#endif
#ifdef _XOPEN_SOURCE
# define __USE_XOPEN 1
# if (_XOPEN_SOURCE - 0) >= 500
# define __USE_XOPEN_EXTENDED 1
# define __USE_UNIX98 1
# undef _LARGEFILE_SOURCE
# define _LARGEFILE_SOURCE 1
# if (_XOPEN_SOURCE - 0) >= 600
# define __USE_XOPEN2K 1
# undef __USE_ISOC99
# define __USE_ISOC99 1
# endif
# else
# ifdef _XOPEN_SOURCE_EXTENDED
# define __USE_XOPEN_EXTENDED 1
# endif
# endif
#endif
#ifdef _LARGEFILE_SOURCE
# define __USE_LARGEFILE 1
#endif
#ifdef _LARGEFILE64_SOURCE
# define __USE_LARGEFILE64 1
#endif
#if defined _FILE_OFFSET_BITS && _FILE_OFFSET_BITS == 64
# define __USE_FILE_OFFSET64 1
#endif
#if defined _BSD_SOURCE || defined _SVID_SOURCE
# define __USE_MISC 1
#endif
#ifdef _BSD_SOURCE
# define __USE_BSD 1
#endif
#ifdef _SVID_SOURCE
# define __USE_SVID 1
#endif
#ifdef _ATFILE_SOURCE
# define __USE_ATFILE 1
#endif
#ifdef _GNU_SOURCE
# define __USE_GNU 1
#endif
#if defined _REENTRANT || defined _THREAD_SAFE
# define __USE_REENTRANT 1
#endif
#if defined _FORTIFY_SOURCE && _FORTIFY_SOURCE > 0 \
&& __GNUC_PREREQ (4, 1) && defined __OPTIMIZE__ && __OPTIMIZE__ > 0
# if _FORTIFY_SOURCE > 1
# define __USE_FORTIFY_LEVEL 2
# else
# define __USE_FORTIFY_LEVEL 1
# endif
#else
# define __USE_FORTIFY_LEVEL 0
#endif
/* We do support the IEC 559 math functionality, real and complex. */
#define __STDC_IEC_559__ 1
#define __STDC_IEC_559_COMPLEX__ 1
/* wchar_t uses ISO 10646-1 (2nd ed., published 2000-09-15) / Unicode 3.1. */
#define __STDC_ISO_10646__ 200009L
/* This macro indicates that the installed library is the GNU C Library.
For historic reasons the value now is 6 and this will stay from now
on. The use of this variable is deprecated. Use __GLIBC__ and
__GLIBC_MINOR__ now (see below) when you want to test for a specific
GNU C library version and use the values in <gnu/lib-names.h> to get
the sonames of the shared libraries. */
#undef __GNU_LIBRARY__
#define __GNU_LIBRARY__ 6
/* Major and minor version number of the GNU C library package. Use
these macros to test for features in specific releases. */
#define __GLIBC__ 2
#define __GLIBC_MINOR__ 4
#define __GLIBC_PREREQ(maj, min) \
((__GLIBC__ << 16) + __GLIBC_MINOR__ >= ((maj) << 16) + (min))
/* Decide whether a compiler supports the long long datatypes. */
#if defined __GNUC__ \
|| (defined __PGI && defined __i386__ ) \
|| (defined __INTEL_COMPILER && (defined __i386__ || defined __ia64__)) \
|| (defined __STDC_VERSION__ && __STDC_VERSION__ >= 199901L)
# define __GLIBC_HAVE_LONG_LONG 1
#endif
/* This is here only because every header file already includes this one. */
#ifndef __ASSEMBLER__
# ifndef _SYS_CDEFS_H
//# include <sys/cdefs.h>
# endif
/* If we don't have __REDIRECT, prototypes will be missing if
__USE_FILE_OFFSET64 but not __USE_LARGEFILE[64]. */
# if defined __USE_FILE_OFFSET64 && !defined __REDIRECT
# define __USE_LARGEFILE 1
# define __USE_LARGEFILE64 1
# endif
#endif /* !ASSEMBLER */
/* Decide whether we can define 'extern inline' functions in headers. */
#if __GNUC_PREREQ (2, 7) && defined __OPTIMIZE__ \
&& !defined __OPTIMIZE_SIZE__ && !defined __NO_INLINE__
# define __USE_EXTERN_INLINES 1
#endif
/* This is here only because every header file already includes this one.
Get the definitions of all the appropriate `__stub_FUNCTION' symbols.
<gnu/stubs.h> contains `#define __stub_FUNCTION' when FUNCTION is a stub
that will always return failure (and set errno to ENOSYS). */
//#include <gnu/stubs.h>
#endif /* features.h */

223
source/shared_lib/include/streflop/libm/headers/ieee754.h Normal file
View File

@@ -0,0 +1,223 @@
/* See the import.pl script for potential modifications */
/* Copyright (C) 1992, 1995, 1996, 1999 Free Software Foundation, Inc.
This file is part of the GNU C Library.
The GNU C Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1 of the License, or (at your option) any later version.
The GNU C Library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with the GNU C Library; if not, write to the Free
Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
02111-1307 USA. */
#ifndef _IEEE754_H
#define _IEEE754_H 1
#include "../streflop_libm_bridge.h"
#include "features.h"
#include "endian.h"
//__BEGIN_DECLS
#if defined(LIBM_COMPILING_FLT32)
union ieee754_float
{
int storage[sizeof(float)/sizeof(int)];
inline Simple& f() {return SIMPLE_FROM_INT_PTR(&storage[0]);}
inline const Simple& f() const {return CONST_SIMPLE_FROM_INT_PTR(&storage[0]);}
/* This is the IEEE 754 single-precision format. */
struct
{
#if __BYTE_ORDER == __BIG_ENDIAN
unsigned int negative:1;
unsigned int exponent:8;
unsigned int mantissa:23;
#endif /* Big endian. */
#if __BYTE_ORDER == __LITTLE_ENDIAN
unsigned int mantissa:23;
unsigned int exponent:8;
unsigned int negative:1;
#endif /* Little endian. */
} ieee;
/* This format makes it easier to see if a NaN is a signalling NaN. */
struct
{
#if __BYTE_ORDER == __BIG_ENDIAN
unsigned int negative:1;
unsigned int exponent:8;
unsigned int quiet_nan:1;
unsigned int mantissa:22;
#endif /* Big endian. */
#if __BYTE_ORDER == __LITTLE_ENDIAN
unsigned int mantissa:22;
unsigned int quiet_nan:1;
unsigned int exponent:8;
unsigned int negative:1;
#endif /* Little endian. */
} ieee_nan;
};
#define IEEE754_FLOAT_BIAS 0x7f /* Added to exponent. */
#endif
#if defined(LIBM_COMPILING_DBL64)
union ieee754_double
{
int storage[sizeof(double)/sizeof(int)];
inline Double& d() {return DOUBLE_FROM_INT_PTR(&storage[0]);}
inline const Double& d() const {return CONST_DOUBLE_FROM_INT_PTR(&storage[0]);}
/* This is the IEEE 754 double-precision format. */
struct
{
#if __BYTE_ORDER == __BIG_ENDIAN
unsigned int negative:1;
unsigned int exponent:11;
/* Together these comprise the mantissa. */
unsigned int mantissa0:20;
unsigned int mantissa1:32;
#endif /* Big endian. */
#if __BYTE_ORDER == __LITTLE_ENDIAN
# if __FLOAT_WORD_ORDER == BIG_ENDIAN
unsigned int mantissa0:20;
unsigned int exponent:11;
unsigned int negative:1;
unsigned int mantissa1:32;
# else
/* Together these comprise the mantissa. */
unsigned int mantissa1:32;
unsigned int mantissa0:20;
unsigned int exponent:11;
unsigned int negative:1;
# endif
#endif /* Little endian. */
} ieee;
/* This format makes it easier to see if a NaN is a signalling NaN. */
struct
{
#if __BYTE_ORDER == __BIG_ENDIAN
unsigned int negative:1;
unsigned int exponent:11;
unsigned int quiet_nan:1;
/* Together these comprise the mantissa. */
unsigned int mantissa0:19;
unsigned int mantissa1:32;
#else
# if __FLOAT_WORD_ORDER == BIG_ENDIAN
unsigned int mantissa0:19;
unsigned int quiet_nan:1;
unsigned int exponent:11;
unsigned int negative:1;
unsigned int mantissa1:32;
# else
/* Together these comprise the mantissa. */
unsigned int mantissa1:32;
unsigned int mantissa0:19;
unsigned int quiet_nan:1;
unsigned int exponent:11;
unsigned int negative:1;
# endif
#endif
} ieee_nan;
};
#define IEEE754_DOUBLE_BIAS 0x3ff /* Added to exponent. */
#endif
#if defined(LIBM_COMPILING_LDBL96)
#if defined(Extended)
union ieee854_long_double
{
int storage[sizeof(long double)/sizeof(int)];
inline Extended& d() {return EXTENDED_FROM_INT_PTR(&storage[0]);}
inline const Extended& d() const {return CONST_EXTENDED_FROM_INT_PTR(&storage[0]);}
/* This is the IEEE 854 double-extended-precision format. */
struct
{
#if __BYTE_ORDER == __BIG_ENDIAN
unsigned int negative:1;
unsigned int exponent:15;
unsigned int empty:16;
unsigned int mantissa0:32;
unsigned int mantissa1:32;
#endif
#if __BYTE_ORDER == __LITTLE_ENDIAN
# if __FLOAT_WORD_ORDER == BIG_ENDIAN
unsigned int exponent:15;
unsigned int negative:1;
unsigned int empty:16;
unsigned int mantissa0:32;
unsigned int mantissa1:32;
# else
unsigned int mantissa1:32;
unsigned int mantissa0:32;
unsigned int exponent:15;
unsigned int negative:1;
unsigned int empty:16;
# endif
#endif
} ieee;
/* This is for NaNs in the IEEE 854 double-extended-precision format. */
struct
{
#if __BYTE_ORDER == __BIG_ENDIAN
unsigned int negative:1;
unsigned int exponent:15;
unsigned int empty:16;
unsigned int one:1;
unsigned int quiet_nan:1;
unsigned int mantissa0:30;
unsigned int mantissa1:32;
#endif
#if __BYTE_ORDER == __LITTLE_ENDIAN
# if __FLOAT_WORD_ORDER == BIG_ENDIAN
unsigned int exponent:15;
unsigned int negative:1;
unsigned int empty:16;
unsigned int mantissa0:30;
unsigned int quiet_nan:1;
unsigned int one:1;
unsigned int mantissa1:32;
# else
unsigned int mantissa1:32;
unsigned int mantissa0:30;
unsigned int quiet_nan:1;
unsigned int one:1;
unsigned int exponent:15;
unsigned int negative:1;
unsigned int empty:16;
# endif
#endif
} ieee_nan;
};
#define IEEE854_LONG_DOUBLE_BIAS 0x3fff
#endif
#endif
//__END_DECLS
#endif /* ieee754.h */

5
source/shared_lib/include/streflop/libm/headers/math.h Normal file
View File

@@ -0,0 +1,5 @@
/* This file is stub automatically generated by import.pl */
// Include bridge, just in case math_private is not included
#include "../streflop_libm_bridge.h"

285
source/shared_lib/include/streflop/libm/headers/math_private.h Normal file
View File

@@ -0,0 +1,285 @@
/* See the import.pl script for potential modifications */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* from: @(#)fdlibm.h 5.1 93/09/24
*/
#ifndef _MATH_PRIVATE_H_
#define _MATH_PRIVATE_H_
/* import.pl: Skipped the macro definitions, keep only the declarations, converted to use streflop types (aliases or wrappers) */
#include "../streflop_libm_bridge.h"
namespace streflop_libm {
#ifdef LIBM_COMPILING_FLT32
#define __sqrtf __ieee754_sqrtf
#define fabsf __fabsf
#define copysignf __copysignf
extern Simple __log1pf(Simple x);
extern Simple __fabsf(Simple x);
extern Simple __atanf(Simple x);
extern Simple __expm1f(Simple x);
extern int __isinff(Simple x);
extern Simple __rintf(Simple x);
extern Simple __cosf(Simple x);
extern void __sincosf (Simple x, Simple *sinx, Simple *cosx);
extern Simple __floorf(Simple x);
extern Simple __scalbnf (Simple x, int n);
extern Simple __frexpf(Simple x, int *eptr);
extern Simple __ldexpf(Simple value, int exp);
extern int __finitef(Simple x);
#endif
#ifdef LIBM_COMPILING_DBL64
#define __sqrt __ieee754_sqrt
#define fabs __fabs
#define copysign __copysign
extern Double __log1p(Double x);
extern Double __fabs(Double x);
extern Double atan(Double x);
extern Double __expm1(Double x);
extern int __isinf(Double x);
extern Double __rint(Double x);
extern Double __cos(Double x);
extern void __sincos (Double x, Double *sinx, Double *cosx);
extern Double __floor(Double x);
extern Double __scalbn(Double x, int n);
extern Double __frexp(Double x, int *eptr);
extern Double __ldexp(Double value, int exp);
extern int __finite(Double x);
#endif
#ifdef LIBM_COMPILING_LDBL96
#if defined(Extended)
#define fabsl __fabsl
extern Extended __cosl(Extended x);
extern Extended __sinl(Extended x);
extern Extended __fabsl(Extended x);
#endif
#endif
#ifdef LIBM_COMPILING_DBL64
extern Double __ieee754_sqrt (Double);
extern Double __ieee754_acos (Double);
extern Double __ieee754_acosh (Double);
extern Double __ieee754_log (Double);
extern Double __ieee754_atanh (Double);
extern Double __ieee754_asin (Double);
extern Double __ieee754_atan2 (Double,Double);
extern Double __ieee754_exp (Double);
extern Double __ieee754_exp2 (Double);
extern Double __ieee754_exp10 (Double);
extern Double __ieee754_cosh (Double);
extern Double __ieee754_fmod (Double,Double);
extern Double __ieee754_pow (Double,Double);
extern Double __ieee754_lgamma_r (Double,int *);
extern Double __ieee754_gamma_r (Double,int *);
extern Double __ieee754_lgamma (Double);
extern Double __ieee754_gamma (Double);
extern Double __ieee754_log10 (Double);
extern Double __ieee754_log2 (Double);
extern Double __ieee754_sinh (Double);
extern Double __ieee754_hypot (Double,Double);
extern Double __ieee754_j0 (Double);
extern Double __ieee754_j1 (Double);
extern Double __ieee754_y0 (Double);
extern Double __ieee754_y1 (Double);
extern Double __ieee754_jn (int,Double);
extern Double __ieee754_yn (int,Double);
extern Double __ieee754_remainder (Double,Double);
extern int32_t __ieee754_rem_pio2 (Double,Double*);
extern Double __ieee754_scalb (Double,Double);
#endif
/* fdlibm kernel function */
#ifdef LIBM_COMPILING_DBL64
extern Double __kernel_standard (Double,Double,int);
extern Double __kernel_sin (Double,Double,int);
extern Double __kernel_cos (Double,Double);
extern Double __kernel_tan (Double,Double,int);
extern int __kernel_rem_pio2 (Double*,Double*,int,int,int, const int32_t*);
#endif
/* internal functions. */
#ifdef LIBM_COMPILING_DBL64
extern Double __copysign (Double x, Double __y);
#endif
#if 0
#ifdef LIBM_COMPILING_DBL64
extern inline Double __copysign (Double x, Double y)
{ return __builtin_copysign (x, y); }
#endif
#endif
/* ieee style elementary Simple functions */
#ifdef LIBM_COMPILING_FLT32
extern Simple __ieee754_sqrtf (Simple);
extern Simple __ieee754_acosf (Simple);
extern Simple __ieee754_acoshf (Simple);
extern Simple __ieee754_logf (Simple);
extern Simple __ieee754_atanhf (Simple);
extern Simple __ieee754_asinf (Simple);
extern Simple __ieee754_atan2f (Simple,Simple);
extern Simple __ieee754_expf (Simple);
extern Simple __ieee754_exp2f (Simple);
extern Simple __ieee754_exp10f (Simple);
extern Simple __ieee754_coshf (Simple);
extern Simple __ieee754_fmodf (Simple,Simple);
extern Simple __ieee754_powf (Simple,Simple);
extern Simple __ieee754_lgammaf_r (Simple,int *);
extern Simple __ieee754_gammaf_r (Simple,int *);
extern Simple __ieee754_lgammaf (Simple);
extern Simple __ieee754_gammaf (Simple);
extern Simple __ieee754_log10f (Simple);
extern Simple __ieee754_log2f (Simple);
extern Simple __ieee754_sinhf (Simple);
extern Simple __ieee754_hypotf (Simple,Simple);
extern Simple __ieee754_j0f (Simple);
extern Simple __ieee754_j1f (Simple);
extern Simple __ieee754_y0f (Simple);
extern Simple __ieee754_y1f (Simple);
extern Simple __ieee754_jnf (int,Simple);
extern Simple __ieee754_ynf (int,Simple);
extern Simple __ieee754_remainderf (Simple,Simple);
extern int32_t __ieee754_rem_pio2f (Simple,Simple*);
extern Simple __ieee754_scalbf (Simple,Simple);
#endif
/* Simple versions of fdlibm kernel functions */
#ifdef LIBM_COMPILING_FLT32
extern Simple __kernel_sinf (Simple,Simple,int);
extern Simple __kernel_cosf (Simple,Simple);
extern Simple __kernel_tanf (Simple,Simple,int);
extern int __kernel_rem_pio2f (Simple*,Simple*,int,int,int, const int32_t*);
#endif
/* internal functions. */
#ifdef LIBM_COMPILING_FLT32
extern Simple __copysignf (Simple x, Simple __y);
#endif
#if 0
#ifdef LIBM_COMPILING_FLT32
extern inline Simple __copysignf (Simple x, Simple y)
{ return __builtin_copysignf (x, y); }
#endif
#endif
#if defined(Extended)
/* ieee style elementary Extended functions */
#ifdef LIBM_COMPILING_LDBL96
extern Extended __ieee754_sqrtl (Extended);
extern Extended __ieee754_acosl (Extended);
extern Extended __ieee754_acoshl (Extended);
extern Extended __ieee754_logl (Extended);
extern Extended __ieee754_atanhl (Extended);
extern Extended __ieee754_asinl (Extended);
extern Extended __ieee754_atan2l (Extended,Extended);
extern Extended __ieee754_expl (Extended);
extern Extended __ieee754_exp2l (Extended);
extern Extended __ieee754_exp10l (Extended);
extern Extended __ieee754_coshl (Extended);
extern Extended __ieee754_fmodl (Extended,Extended);
extern Extended __ieee754_powl (Extended,Extended);
extern Extended __ieee754_lgammal_r (Extended,int *);
extern Extended __ieee754_gammal_r (Extended,int *);
extern Extended __ieee754_lgammal (Extended);
extern Extended __ieee754_gammal (Extended);
extern Extended __ieee754_log10l (Extended);
extern Extended __ieee754_log2l (Extended);
extern Extended __ieee754_sinhl (Extended);
extern Extended __ieee754_hypotl (Extended,Extended);
extern Extended __ieee754_j0l (Extended);
extern Extended __ieee754_j1l (Extended);
extern Extended __ieee754_y0l (Extended);
extern Extended __ieee754_y1l (Extended);
extern Extended __ieee754_jnl (int,Extended);
extern Extended __ieee754_ynl (int,Extended);
extern Extended __ieee754_remainderl (Extended,Extended);
extern int __ieee754_rem_pio2l (Extended,Extended*);
extern Extended __ieee754_scalbl (Extended,Extended);
#endif
/* Extended versions of fdlibm kernel functions */
#ifdef LIBM_COMPILING_LDBL96
extern Extended __kernel_sinl (Extended,Extended,int);
extern Extended __kernel_cosl (Extended,Extended);
extern Extended __kernel_tanl (Extended,Extended,int);
extern void __kernel_sincosl (Extended,Extended,
Extended *,Extended *, int);
extern int __kernel_rem_pio2l (Extended*,Extended*,int,int,
int,const int*);
#endif
#ifndef NO_LONG_DOUBLE
/* prototypes required to compile the ldbl-96 support without warnings */
#ifdef LIBM_COMPILING_LDBL96
extern int __finitel (Extended);
extern int __ilogbl (Extended);
extern int __isinfl (Extended);
extern int __isnanl (Extended);
extern Extended __atanl (Extended);
extern Extended __copysignl (Extended, Extended);
extern Extended __expm1l (Extended);
extern Extended __floorl (Extended);
extern Extended __frexpl (Extended, int *);
extern Extended __ldexpl (Extended, int);
extern Extended __log1pl (Extended);
extern Extended __nanl (const char *);
extern Extended __rintl (Extended);
extern Extended __scalbnl (Extended, int);
extern Extended __sqrtl (Extended x);
extern Extended fabsl (Extended x);
extern void __sincosl (Extended, Extended *, Extended *);
extern Extended __logbl (Extended x);
extern Extended __significandl (Extended x);
#endif
#if 0
#ifdef LIBM_COMPILING_LDBL96
extern inline Extended __copysignl (Extended x, Extended y)
{ return __builtin_copysignl (x, y); }
#endif
#endif
#endif
#endif
/* Prototypes for functions of the IBM Accurate Mathematical Library. */
#ifdef LIBM_COMPILING_DBL64
extern Double __exp1 (Double __x, Double __xx, Double __error);
extern Double __sin (Double __x);
extern Double __cos (Double __x);
extern int __branred (Double __x, Double *__a, Double *__aa);
extern void __doasin (Double __x, Double __dx, Double __v[]);
extern void __dubsin (Double __x, Double __dx, Double __v[]);
extern void __dubcos (Double __x, Double __dx, Double __v[]);
extern Double __halfulp (Double __x, Double __y);
extern Double __sin32 (Double __x, Double __res, Double __res1);
extern Double __cos32 (Double __x, Double __res, Double __res1);
extern Double __mpsin (Double __x, Double __dx);
extern Double __mpcos (Double __x, Double __dx);
extern Double __mpsin1 (Double __x);
extern Double __mpcos1 (Double __x);
extern Double __slowexp (Double __x);
extern Double __slowpow (Double __x, Double __y, Double __z);
extern void __docos (Double __x, Double __dx, Double __v[]);
#endif
}
#endif /* _MATH_PRIVATE_H_ */

26
source/shared_lib/include/streflop/libm/headers/wchar.h Normal file
View File

@@ -0,0 +1,26 @@
/* wchar_t type related definitions.
Copyright (C) 2000 Free Software Foundation, Inc.
This file is part of the GNU C Library.
The GNU C Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1 of the License, or (at your option) any later version.
The GNU C Library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with the GNU C Library; if not, write to the Free
Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
02111-1307 USA. */
#ifndef _BITS_WCHAR_H
#define _BITS_WCHAR_H 1
#define __WCHAR_MIN (-2147483647 - 1)
#define __WCHAR_MAX (2147483647)
#endif /* bits/wchar.h */

536
source/shared_lib/include/streflop/libm/import.pl Executable file
View File

@@ -0,0 +1,536 @@
#!/usr/bin/perl
# streflop: STandalone REproducible FLOating-Point
# This script imports the GNU libm files and convert them in such a way they can be compiled with streflop types.
# Code released according to the GNU Lesser General Public License
# Nicolas Brodu, 2006
# Please read the history and copyright information in the accompanying documentation
if ($#ARGV==-1) {
print <<ENDHELP;
Perl script to import the libm components from the GNU glibc inside this project.
Tested with glibc-2.4
Usage: ./import.pl glibc_src_dir
Where: glibc_src_dir is the source directory obtained by uncompressing the GNU libc archive.
This script copies the libm directories for generic IEEE754 math with float (32 bits), double (64 bits) and long double (80 bits). The files are then post-processed to remove all dependencies from glibc so they can be compiled standalone.
ENDHELP
exit 0;
}
$dir=$ARGV[0];
if ((! -d "$dir/sysdeps/ieee754/flt-32")
|| (! -d "$dir/sysdeps/ieee754/dbl-64")
|| (! -d "$dir/sysdeps/ieee754/ldbl-96")
|| (! -r "$dir/math/s_ldexpf.c") # was put here for unknown reason, probably because not system-dependent
|| (! -r "$dir/math/s_ldexp.c") # was put here for unknown reason, probably because not system-dependent
|| (! -r "$dir/math/s_nextafter.c") # was put here for unknown reason, probably because not system-dependent
|| (! -r "$dir/math/math_private.h")
|| (! -r "$dir/sysdeps/ieee754/ieee754.h")
|| (! -r "$dir/include/features.h")
|| (! -r "$dir/string/endian.h")
|| (! -r "$dir/include/endian.h")
) {
print "Invalid glibc directory\n";
exit 1;
}
$importNotice="/* See the import.pl script for potential modifications */\n";
# Copy files
system("cp -r $dir/sysdeps/ieee754/flt-32 $dir/sysdeps/ieee754/dbl-64 $dir/sysdeps/ieee754/ldbl-96 .");
if (! -d "headers") {mkdir("headers");}
system("cp -r $dir/math/s_ldexpf.c flt-32");
system("cp -r $dir/math/s_ldexp.c $dir/math/s_nextafter.c dbl-64");
# this ldbl function is actually an alias to the dbl one in the libm. Make a hard copy here.
system("cp -r $dir/math/s_ldexp.c ldbl-96/s_ldexpl.c");
system("cp -r $dir/math/math_private.h $dir/sysdeps/ieee754/ieee754.h $dir/include/features.h $dir/string/endian.h $dir/bits/wchar.h headers/");
# Merge include files
system("cat $dir/include/endian.h >> headers/endian.h");
# remove some routines we don't need
if (-r "ldbl-96/printf_fphex.c") {unlink "ldbl-96/printf_fphex.c";}
if (-r "ldbl-96/strtold_l.c") {unlink "ldbl-96/strtold_l.c";}
if (-r "flt-32/mpn2flt.c") {unlink "flt-32/mpn2flt.c";}
if (-r "dbl-64/mpn2dbl.c") {unlink "dbl-64/mpn2dbl.c";}
if (-r "dbl-64/dbl2mpn.c") {unlink "dbl-64/dbl2mpn.c";}
if (-r "dbl-64/t_exp.c") {unlink "dbl-64/t_exp.c";}
if (-r "ldbl-96/mpn2ldbl.c") {unlink "ldbl-96/mpn2ldbl.c";}
if (-r "ldbl-96/ldbl2mpn.c") {unlink "ldbl-96/ldbl2mpn.c";}
if (-r "ldbl-96/s_nexttoward.c") {unlink "ldbl-96/s_nexttoward.c";}
if (-r "ldbl-96/s_nexttowardf.c") {unlink "ldbl-96/s_nexttowardf.c";}
if (-r "ldbl-96/math_ldbl.h") {unlink "ldbl-96/math_ldbl.h";}
# The float exp in e_expf.c uses doubles internally since revision 1.2 in the libm-ieee754 CVS attic!
# Roll back to the slower, but purely float version, and also overwrite the wrapper by a wraper to the float only version
system("cp -f w_expf.c e_expf.c flt-32");
# convert .c => .cpp for clarity
@filelist = glob("flt-32/*.c dbl-64/*.c ldbl-96/*.c");
foreach $f (@filelist) {
$g = $f;
$g =~ s/\.c$/.cpp/;
rename $f, $g;
}
# create dummy math.h stub
open(FILE, ">headers/math.h");
print FILE <<MATH_H;
/* This file is stub automatically generated by import.pl */
// Include bridge, just in case math_private is not included
#include "../streflop_libm_bridge.h"
MATH_H
close FILE;
# These files in the flt directory use double but could use floats instead
foreach $f ("flt-32/s_cbrtf.cpp", "flt-32/s_llrintf.cpp", "flt-32/s_lrintf.cpp") {
open(FILE,"<$f");
$content = "";
while(<FILE>) {
# replace double by float, OK in these files
s/double/float/g;
$content.=$_;
}
close FILE;
open(FILE,">$f");
print FILE $content;
close FILE;
}
# These files in the dbl directory wrongly use float, long double and float functions
foreach $f ("dbl-64/e_lgamma_r.cpp", "dbl-64/t_exp2.h", "dbl-64/s_llrint.cpp") {
open(FILE,"<$f");
$content = "";
while(<FILE>) {
s/fabsf/fabs/g;
s/float/double/g;
s/long double/double/g;
$content.=$_;
}
close FILE;
open(FILE,">$f");
print FILE $content;
close FILE;
}
# These files in the ldbl directory wrongly use lower precision types
foreach $f ("ldbl-96/e_lgammal_r.cpp", "ldbl-96/s_cbrtl.cpp", "ldbl-96/s_logbl.cpp", "ldbl-96/s_ldexpl.cpp") {
open(FILE,"<$f");
$content = "";
while(<FILE>) {
# do the substitution now
s/long double/Extended/g;
s/\(double\) i;/(long double) i;/g;
s/double (_|v)/long double $1/g;
s/const double factor/const long double factor/g;
s/fabs(?!l)/fabsl$1/g;
s/ldexp(?!l)/ldexpl$1/g;
s/__finite\(/__finitel\(/g;
s/__scalbn\(/__scalbnl\(/g;
$content.=$_;
}
close FILE;
open(FILE,">$f");
print FILE $content;
close FILE;
}
# DOUBLE_FROM_INT_PTR(x) could simply be *reinterpret_cast<double*>(x)
# for basic types, but may use a dedicated factory for Object wrappers
# BaseType is either the same as FloatType for plain old float/double, or it is
# the float/double type when FloatType is an Object wrapper
$xdaccessor=
"inline Double& d() {return DOUBLE_FROM_INT_PTR(&i[0]);}\n"
."inline Double& x() {return DOUBLE_FROM_INT_PTR(&i[0]);}\n"
."inline Double& d(int idx) {return DOUBLE_FROM_INT_PTR(&i[idx*(sizeof(double)/sizeof(i))]);}\n"
."inline Double& x(int idx) {return DOUBLE_FROM_INT_PTR(&i[idx*(sizeof(double)/sizeof(i))]);}\n"
."inline const Double& d() const {return CONST_DOUBLE_FROM_INT_PTR(&i[0]);}\n"
."inline const Double& x() const {return CONST_DOUBLE_FROM_INT_PTR(&i[0]);}\n"
."inline const Double& d(int idx) const {return CONST_DOUBLE_FROM_INT_PTR(&i[idx*(sizeof(double)/sizeof(i))]);}\n"
."inline const Double& x(int idx) const {return CONST_DOUBLE_FROM_INT_PTR(&i[idx*(sizeof(double)/sizeof(i))]);}\n"
;
@filelist = glob("flt-32/* dbl-64/* ldbl-96/*");
foreach $f (@filelist) {
open(FILE,"<$f");
$content = "";
$opened_namespace = 0;
while(<FILE>) {
# remove all system-wide includes
if (/.*?#include(.*?)<sys\/(.*)/) {$_="//".$_;}
# Convert all includes to local files
s/^([ \t]*)#include([ \t]*)<(.*)>(.*)/#$1include$2\"$3\"$4/g;
# fp environment is now handled directly by streflop (and FPUSettings.h)
s/fenv\.h/..\/streflop_libm_bridge.h/g;
# integer types handled by the bridge
s/inttypes\.h/..\/streflop_libm_bridge.h/g;
# get rid of stdlib too
s/stdlib\.h/..\/streflop_libm_bridge.h/g;
# and limits.h
s/limits\.h/..\/streflop_libm_bridge.h/g;
# float.h too
s/float\.h/..\/streflop_libm_bridge.h/g;
# errno.h falls back to system, remove it
s/errno\.h/..\/streflop_libm_bridge.h/g;
# replace all occurences of problematic union fields with objects
# by proper C++ accessors
s/\b(\.d(?!\[)\b|\.d\[(.*?)\])/.d($2)/g;
s/\b(\.x(?!\[)\b|\.x\[(.*?)\])/.x($2)/g;
s/\b(\.f(?!\[)\b|\.f\[(.*?)\])/.f($2)/g; # ieee754.h
s/\b(->d(?!\[)\b|->d\[(.*?)\])/->d($2)/g;
s/\b(->x(?!\[)\b|->x\[(.*?)\])/->x($2)/g;
s/\b(->f(?!\[)\b|->f\[(.*?)\])/->f($2)/g; # ieee754.h
# Some more occurences from arrays of such unions/structs
s/\]\.d\b/].d()/g;
# named field C construct is invalid (C++ would initialize first union member)
# and we replace unions by struct + accessor anyway
s/{ *?.i *?= *?{/{{/g;
# volatile declaration of static const global variables poses problem with the wrapper types
s/volatile //g;
# Now substitute the base types by their Simple/Double/Extended aliases or wrapper
s/long double/Extended/g; # before double
s/\bdouble\b/Double/g;
s/\bfloat\b/Simple/g;
# Replace problematic int += double
s/E(X|Y|Z)(.*?=.*?)ONE/E${1}${2}1/g;
# problematic ?: operator with different types. This simple check catches all problematic occurences
# if (/\?(.*?(\b0\.|\.0).*?:.*?|.*?:.*?(\b0\.|\.0).*?);/) {
# print "$f => ?$1;\n";
# }
if (/\?(.*?):(.*?)\b(0\.0|10\.0|1\.0)\b/) {
my $type = "";
if ($f =~ /flt-32/) {$type = "Simple";}
if ($f =~ /dbl-64/) {$type = "Double";}
if ($f =~ /ldbl-96/) {$type = "Extended";}
s/\?(.*?):(.*?)\b(0\.0|10\.0|1\.0)\b/?$1:$2$type($3)/g;
}
my $type = "";
my $flit = "";
if ($f =~ /flt-32/) {$type = "Simple"; $flit = "f";}
if ($f =~ /dbl-64/) {$type = "Double"; $flit = "";}
if ($f =~ /ldbl-96/) {$type = "Extended"; $flit = "l";}
# replace problematic mixed-type ?: operators where an int and a float are used as arguments, incompatible with wrappers
if (/\?(.*?)\b(0\.0|1\.0)\b(.*?):(.*?)/) {
s/\?(.*?)\b(0\.0|1\.0)\b(.*?):(.*?)/?$1$type($2)$3:$4/g;
}
# These special cases are OK because no other ?: pattern use them in the 3 subdirs
s/\?0:/?Double(0.0):/;
s/:0;/:Double(0.0);/;
# protect the new symbol names by namespace to avoid any conflict with system libm
if (((/#ifdef __STDC__/) || (/.*? (__|mcr|ss32)[a-z,A-Z,_,0-9]*?\(.*?{$/) || (/Double (atan2Mp|atanMp|slow|tanMp|__exp1|__ieee754_remainder|__ieee754_sqrt)/) || (/^(Simple|Double|Extended|void|int)$/) || ((/^#ifdef BIG_ENDI$/) && ($f =~ /(uatan|mpa2|mpexp|atnat|sincos32)/)) || (/^#define MM 5$/) || (/^void __mp(log|sqrt|exp|atan)\(/) || (/^int __(b|mp)ranred\(/)) && ($opened_namespace == 0)) {
$_ = "namespace streflop_libm {\n".$_;
$opened_namespace = 1;
}
# Prevent type promotion for native aliases
# Ex: promotion would cause computations like y = 0.5 * x to be done in double for x,y float, then cast back to float
# We want the whole computation done in float, not temporary switching to double.
# In some case the FPU internal precision is enough to mask the problem, but for SSE for example, the float/double
# are using different instructions and the problem cannot be ignored (in these cases, there are differences from
# soft float implementation)
# Cannot replace by s/blah/$type($1)/g; because of static order initialization fiasco with wrapper types
# So use next best option to use litteral type specification
# => this solves the problem for native types
# => wrappers are OK thanks to the redefinitions of the operators
s/\b((\d+\.\d*|\d*\.\d+)((e|E)[-+]?\d+)?)(f|F|l|L)?\b/$1$flit/g;
$content.=$_;
}
close FILE;
# multi-line spanning regexp
$content =~ s/union ?{(.*?);.*?Double.*?}/struct {\n$xdaccessor$1;}/sg;
# close opened namespace
if ($opened_namespace==1) {
$content .= "}\n";
}
open(FILE,">$f");
print FILE $importNotice;
print FILE $content;
close FILE;
}
# ieee754.h union+accessor
open(FILE,"<headers/ieee754.h");
$content = "";
while(<FILE>) {
# Comment out decls sections
if (/.*?__.*?_DECLS.*/) {$_="//".$_;}
# Convert all includes to local files
s/^([ \t]*)#include([ \t]*)<(.*)>(.*)/#$1include$2\"$3\"$4/g;
# insert our own bridge
if (/^#define _IEEE754_H.*/) {
$_.="#include \"../streflop_libm_bridge.h\"\n\n";
}
# Protect the Simple section by a #define
if (/.*?ieee754_float.*?/) {
$_ = "#if defined(LIBM_COMPILING_FLT32)\n".$_;
}
if (/.*?IEEE754_FLOAT_BIAS.*?/) {
$_ = $_."\n#endif\n";
}
# Protect the Double section by a #define
if (/.*?ieee754_double.*?/) {
$_ = "#if defined(LIBM_COMPILING_DBL64)\n".$_;
}
if (/.*?IEEE754_DOUBLE_BIAS.*?/) {
$_ = $_."\n#endif\n";
}
# Protect the Extended section by a #define
if (/.*?ieee854_long_double.*?/) {
$_ = "#if defined(Extended)\n".$_;
$_ = "#if defined(LIBM_COMPILING_LDBL96)\n".$_;
}
if (/.*?IEEE854_LONG_DOUBLE_BIAS.*?/) {
$_ = $_."\n#endif\n";
$_ = $_."\n#endif\n";
}
$content.=$_;
}
$ieeeAccessorSimple=
"inline Simple& f() {return SIMPLE_FROM_INT_PTR(&storage[0]);}\n"
."inline const Simple& f() const {return CONST_SIMPLE_FROM_INT_PTR(&storage[0]);}\n";
$ieeeAccessorDouble=
"inline Double& d() {return DOUBLE_FROM_INT_PTR(&storage[0]);}\n"
."inline const Double& d() const {return CONST_DOUBLE_FROM_INT_PTR(&storage[0]);}\n";
$ieeeAccessorExtended=
"inline Extended& d() {return EXTENDED_FROM_INT_PTR(&storage[0]);}\n"
."inline const Extended& d() const {return CONST_EXTENDED_FROM_INT_PTR(&storage[0]);}\n";
# multi-line spanning regexp
$content =~ s/union(.*?ieee854_long_double.*?){.*?;/union$1\{\nint storage[sizeof(long double)\/sizeof(int)];\n$ieeeAccessorExtended/sg;
$content =~ s/union(.*?ieee754_double.*?){.*?;/union$1\{\nint storage[sizeof(double)\/sizeof(int)];\n$ieeeAccessorDouble/sg;
$content =~ s/union(.*?ieee754_float.*?){.*?;/union$1\{\nint storage[sizeof(float)\/sizeof(int)];\n$ieeeAccessorSimple/sg;
close FILE;
open(FILE,">headers/ieee754.h");
print FILE $importNotice;
print FILE $content;
close FILE;
# math_private.h needs a special treatement to define the macros for the wrapper objects
# The macros will be defined separately
open(FILE,"<headers/math_private.h");
my $flag=1;
my $precisionMode = "none";
@convert=();
MPRIV_LOOP: while(<FILE>) {
# keep initial lines with header information and #define file protection
# skip all lines defining macros and unions
if (/^#define _MATH_PRIVATE_H_.*/) {
push @convert,$_;
push @convert,"\n/* import.pl: Skipped the macro definitions, keep only the declarations, converted to use streflop types (aliases or wrappers) */\n#include \"../streflop_libm_bridge.h\"\n\nnamespace streflop_libm {\n\n";
# ensure strict separation of precision mode, each one does not have access to the other names
push @convert,"#ifdef LIBM_COMPILING_FLT32\n";
# define wrappers as the libm would do in IEEE754 mode
push @convert,"#define __sqrtf __ieee754_sqrtf\n";
push @convert,"#define fabsf __fabsf\n";
push @convert,"#define copysignf __copysignf\n";
# add missing defines
push @convert,"extern Simple __log1pf(Simple x);\n";
push @convert,"extern Simple __fabsf(Simple x);\n";
push @convert,"extern Simple __atanf(Simple x);\n";
push @convert,"extern Simple __expm1f(Simple x);\n";
push @convert,"extern int __isinff(Simple x);\n";
push @convert,"extern Simple __rintf(Simple x);\n";
push @convert,"extern Simple __cosf(Simple x);\n";
push @convert,"extern void __sincosf (Simple x, Simple *sinx, Simple *cosx);\n";
push @convert,"extern Simple __floorf(Simple x);\n";
push @convert,"extern Simple __scalbnf (Simple x, int n);\n";
push @convert,"extern Simple __frexpf(Simple x, int *eptr);\n";
push @convert,"extern Simple __ldexpf(Simple value, int exp);\n";
push @convert,"extern int __finitef(Simple x);\n";
push @convert,"#endif\n\n";
push @convert,"#ifdef LIBM_COMPILING_DBL64\n";
push @convert,"#define __sqrt __ieee754_sqrt\n";
push @convert,"#define fabs __fabs\n";
push @convert,"#define copysign __copysign\n";
push @convert,"extern Double __log1p(Double x);\n";
push @convert,"extern Double __fabs(Double x);\n";
push @convert,"extern Double atan(Double x);\n";
push @convert,"extern Double __expm1(Double x);\n";
push @convert,"extern int __isinf(Double x);\n";
push @convert,"extern Double __rint(Double x);\n";
push @convert,"extern Double __cos(Double x);\n";
push @convert,"extern void __sincos (Double x, Double *sinx, Double *cosx);\n";
push @convert,"extern Double __floor(Double x);\n";
push @convert,"extern Double __scalbn(Double x, int n);\n";
push @convert,"extern Double __frexp(Double x, int *eptr);\n";
push @convert,"extern Double __ldexp(Double value, int exp);\n";
push @convert,"extern int __finite(Double x);\n";
push @convert,"#endif\n\n";
push @convert,"#ifdef LIBM_COMPILING_LDBL96\n";
push @convert,"#if defined(Extended)\n";
push @convert,"#define fabsl __fabsl\n";
push @convert,"extern Extended __cosl(Extended x);\n";
push @convert,"extern Extended __sinl(Extended x);\n";
push @convert,"extern Extended __fabsl(Extended x);\n";
push @convert,"#endif\n";
push @convert,"#endif\n";
push @convert,"\n";
$flag = 0;
}
if (/^extern(.*)/) {
$flag = 1;
}
if ($flag==0) {next MPRIV_LOOP;}
# Now substitute the base types by their Simple/Double/Extended aliases or wrapper
s/\blong double\b/Extended/g; # before double
s/\bdouble\b/Double/g;
s/\bfloat\b/Simple/g;
# Protect the Extended section by a #define
if (/.*?elementary Extended functions.*?/) {
$_ = "#if defined(Extended)\n".$_;
}
if (/.*?functions of the IBM.*?/) {
$_ = "#endif\n\n".$_;
}
# remove the inline aliases
s/__GNUC_PREREQ \(.*?\)/0/;
# end namespace protection
if (/^#endif.*_MATH_PRIVATE_H_/) {
$_ = "}\n\n".$_;
}
# now insert protection for symbols separation
if (/^extern(.*)/) {
$remline = $1;
if ($remline =~ /Extended/) {
if ($precisionMode ne "Extended") {
if ($precisionMode ne "none") {$_ = "#endif\n".$_;}
$_ = "#ifdef LIBM_COMPILING_LDBL96\n".$_;
$precisionMode = "Extended";
}
}
elsif ($remline =~ /Double/) {
if ($precisionMode ne "Double") {
if ($precisionMode ne "none") {$_ = "#endif\n".$_;}
$_ = "#ifdef LIBM_COMPILING_DBL64\n".$_;
$precisionMode = "Double";
}
}
elsif ($remline =~ /Simple/) {
if ($precisionMode ne "Simple") {
if ($precisionMode ne "none") {$_ = "#endif\n".$_;}
$_ = "#ifdef LIBM_COMPILING_FLT32\n".$_;
$precisionMode = "Simple";
}
}
} else {
$line = $_;
chomp $line;
if (($line =~ /^(\s)*$/) && ($precisionMode ne "none")) {
$_ = "#endif\n".$_;
$precisionMode = "none";
}
}
push @convert,$_;
}
close FILE;
open(FILE,">headers/math_private.h");
print FILE $importNotice;
print FILE @convert;
close FILE;
# features.h is nearly ready, just do not include more external defs
open(FILE,"<headers/features.h");
@convert=();
while(<FILE>) {
# commment out external includes
s,(.*?#.*?include.*),//$1,;
push @convert,$_;
}
close FILE;
open(FILE,">headers/features.h");
print FILE $importNotice;
print FILE @convert;
close FILE;
# headers/endian.h
open(FILE,"<headers/endian.h");
@convert=();
while(<FILE>) {
# change machine specific bits/endian by streflop configuration
if (/<bits\/endian\.h>/) {
s,<bits/endian\.h>,\"../streflop_libm_bridge.h\",;
#keep features.h, but locally
} elsif (/<features\.h>/) {
s,<features\.h>,"features.h",;
} else {
# commment out all other external includes
s,(.*?#.*?include.*),//$1,;
}
#unconditional definition of endian things
if (/defined _LIBC/) {
$_ = "#if 1\n//".$_;
}
push @convert,$_;
}
close FILE;
open(FILE,">headers/endian.h");
print FILE $importNotice;
print FILE @convert;
close FILE;
# include the bridge from mpa.h
open(FILE,"<dbl-64/mpa.h");
$content="#include \"../streflop_libm_bridge.h\"\nnamespace streflop_libm {";
while(<FILE>) {
# special case
s/->d\(\)/->mantissa/g;
$content.=$_;
}
close FILE;
$mpAccessor=
"inline Double& d(int idx) {return mantissa[idx];}\n"
."inline const Double& d(int idx) const {return mantissa[idx];}\n";
# multi-line spanning regexp
$content =~ s/Double d\[(.*?)\].*?} mp_no/Double mantissa[$1];\n$mpAccessor} mp_no/sg;
open(FILE,">dbl-64/mpa.h");
print FILE $content."}\n"; # also close namespace
close FILE;
# Generate specific Makefiles
$parentMakefile="";
foreach $dir ("flt-32", "dbl-64", "ldbl-96") {
@objFiles = glob("$dir/*.cpp");
$parentObjects = "$dir-objects =";
foreach $f (@objFiles) {
$f =~ s/\.cpp$/.o/;
$parentObjects .= " libm/$f";
$f =~ s/^$dir\///;
}
open(FILE, ">$dir/Makefile");
my $DIR = uc($dir);
$DIR =~ s/\-//;
print FILE "# Makefile automatically generated by import.pl\n"
."include ../../Makefile.common\n"
."CPPFLAGS += -I../headers -DLIBM_COMPILING_$DIR=1\n"
."all: @objFiles\n"
."\techo '$dir done!'\n";
close FILE;
$parentMakefile .= $parentObjects."\n\n";
}
# generate Makefile rule in parent directory
@projectFiles = glob("flt-32/* dbl-64/* ldbl-96/* headers/*");
foreach $f (@projectFiles) {$f = " libm/$f";}
open(FILE, ">../Makefile.libm_objects");
print FILE "# Makefile automatically generated by libm/import.pl. Do not edit.\n\n"
.$parentMakefile
."\nlibm-src =";
print FILE @projectFiles;
print FILE "\n";
close FILE;

438
source/shared_lib/include/streflop/libm/streflop_libm_bridge.h Normal file
View File

@@ -0,0 +1,438 @@
/*
This file acts as a bridge between the libm files and streflop.
It replaces the missing type definitions and macros that were defined
by the external files that were included by the libm, so it is now
standalone.
In addition, this allows to wrap the basic types for compiling
the libm, ex, with SoftFloat.
This file is NOT included by the streflop main files, and thus NOT
included by the user programs.
It is included by the libm converted files to produce the streflop binary,
that will be linked by the user program instead of the system libm.
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
*/
#ifndef STREFLOP_LIBM_BRIDGE
#define STREFLOP_LIBM_BRIDGE
// prevent inclusion of the main Math.h file
// => the function symbol are not defined when compiling libm
// => this adds a level of protection, a function inadvertantly using a wrong precision function is detected
#define STREFLOP_MATH_H
// First define our custom types
#include "../streflop.h"
// Export type definitions, etc.
using namespace streflop;
// for int32_t, etc. Assume WORDSIZE is defined by the Makefile.common
#include "../System.h"
#undef __USE_ISOC99
#undef __USE_POSIX
#undef __USE_POSIX2
#undef __USE_POSIX199309
#undef __USE_POSIX199506
#undef __USE_XOPEN
#undef __USE_XOPEN_EXTENDED
#undef __USE_UNIX98
#undef __USE_XOPEN2K
#undef __USE_LARGEFILE
#undef __USE_LARGEFILE64
#undef __USE_FILE_OFFSET64
#undef __USE_BSD
#undef __USE_SVID
#undef __USE_MISC
#undef __USE_GNU
#undef __USE_REENTRANT
#undef __USE_FORTIFY_LEVEL
#undef __FAVOR_BSD
#undef __KERNEL_STRICT_NAMES
#undef weak_alias
#undef strong_alias
#undef hidden_def
#undef libm_hidden_def
#undef INTDEF
#define weak_alias(x,y)
#define strong_alias(x,y)
#define hidden_def(x)
#define libm_hidden_def(x)
#define INTDEF(x)
#define __set_errno(x) do {} while (0)
#define _IEEE_LIBM 1
// u_int32_t is not C99 compliant
typedef streflop::uint32_t u_int32_t;
/////////////////////////////////////////////////////////////////////////////
// Common definitions
/////////////////////////////////////////////////////////////////////////////
// <float.h> is not included, define only these constants that are needed to compile libm
#define FLT_MIN_EXP (-128)
#define FLT_MAX_EXP 127
#define FLT_MANT_DIG 24
#define DBL_MIN_EXP (-1024)
#define DBL_MAX_EXP 1023
#define DBL_MANT_DIG 53
#define LDBL_MANT_DIG 64
#define HUGE_VALF SimplePositiveInfinity
#define HUGE_VAL DoublePositiveInfinity
#define HUGE_VALL ExtendedPositiveInfinity
#define INT_MAX ((-1)>>1)
// Stolen from bits/mathdef.h
#define FP_ILOGB0 (-2147483647)
#define FP_ILOGBNAN 2147483647
namespace streflop_libm {
// No template, only allow this for int types (could use traits, but it's simpler this way)
inline char abs(char x) {return x < 0 ? -x : x;}
inline short abs(short x) {return x < 0 ? -x : x;}
inline int abs(int x) {return x < 0 ? -x : x;}
inline long abs(long x) {return x < 0 ? -x : x;}
inline long long abs(long long x) {return x < 0 ? -x : x;}
inline char MIN(char x, char y) {return x<y ? x : y;}
inline short MIN(short x, short y) {return x<y ? x : y;}
inline int MIN(int x, int y) {return x<y ? x : y;}
inline long MIN(long x, long y) {return x<y ? x : y;}
inline long long MIN(long long x, long long y) {return x<y ? x : y;}
// Compensate lack of Math.h
enum
{
FP_NAN = 0,
FP_INFINITE = 1,
FP_ZERO = 2,
FP_SUBNORMAL = 3,
FP_NORMAL = 4
};
#ifdef LIBM_COMPILING_FLT32
extern int __fpclassifyf(Simple x);
extern int __isnanf(Simple x);
extern int __isinff(Simple x);
inline int fpclassify(Simple x) {return streflop_libm::__fpclassifyf(x);}
inline int isnan(Simple x) {return streflop_libm::__isnanf(x);}
inline int isinf(Simple x) {return streflop_libm::__isinff(x);}
inline int isfinite(Simple x) {return !(isnan(x) || isinf(x));}
inline int isnormal(Simple x) {return fpclassify(x) == FP_NORMAL;}
inline bool isunordered(Simple x, Simple y) {
return (fpclassify(x) == FP_NAN) || (fpclassify (y) == FP_NAN);
}
inline bool isgreater(Simple x, Simple y) {
return (!isunordered(x,y)) && (x > y);
}
inline bool isgreaterequal(Simple x, Simple y) {
return (!isunordered(x,y)) && (x >= y);
}
inline bool isless(Simple x, Simple y) {
return (!isunordered(x,y)) && (x < y);
}
inline bool islessequal(Simple x, Simple y) {
return (!isunordered(x,y)) && (x <= y);
}
inline bool islessgreater(Simple x, Simple y) {
return (!isunordered(x,y)) && ((x < y) || (x > y));
}
extern const Simple SimplePositiveInfinity;
extern const Simple SimpleNegativeInfinity;
extern const Simple SimpleNaN;
#endif
#ifdef LIBM_COMPILING_DBL64
extern int __fpclassify(Double x);
extern int __isnanl(Double x);
extern int __isinf(Double x);
inline int fpclassify(Double x) {return streflop_libm::__fpclassify(x);}
inline int isnan(Double x) {return streflop_libm::__isnanl(x);}
inline int isinf(Double x) {return streflop_libm::__isinf(x);}
inline int isfinite(Double x) {return !(isnan(x) || isinf(x));}
inline int isnormal(Double x) {return fpclassify(x) == FP_NORMAL;}
inline bool isunordered(Double x, Double y) {
return (fpclassify(x) == FP_NAN) || (fpclassify (y) == FP_NAN);
}
inline bool isgreater(Double x, Double y) {
return (!isunordered(x,y)) && (x > y);
}
inline bool isgreaterequal(Double x, Double y) {
return (!isunordered(x,y)) && (x >= y);
}
inline bool isless(Double x, Double y) {
return (!isunordered(x,y)) && (x < y);
}
inline bool islessequal(Double x, Double y) {
return (!isunordered(x,y)) && (x <= y);
}
inline bool islessgreater(Double x, Double y) {
return (!isunordered(x,y)) && ((x < y) || (x > y));
}
extern const Double DoublePositiveInfinity;
extern const Double DoubleNegativeInfinity;
extern const Double DoubleNaN;
#endif
#ifdef LIBM_COMPILING_LDBL96
#ifdef Extended
extern int __fpclassifyl(Extended x);
extern int __isnanl(Extended x);
extern int __isinfl(Extended x);
inline int fpclassify(Extended x) {return streflop_libm::__fpclassifyl(x);}
inline int isnan(Extended x) {return streflop_libm::__isnanl(x);}
inline int isinf(Extended x) {return streflop_libm::__isinfl(x);}
inline int isfinite(Extended x) {return !(isnan(x) || isinf(x));}
inline int isnormal(Extended x) {return fpclassify(x) == FP_NORMAL;}
inline bool isunordered(Extended x, Extended y) {
return (fpclassify(x) == FP_NAN) || (fpclassify (y) == FP_NAN);
}
inline bool isgreater(Extended x, Extended y) {
return (!isunordered(x,y)) && (x > y);
}
inline bool isgreaterequal(Extended x, Extended y) {
return (!isunordered(x,y)) && (x >= y);
}
inline bool isless(Extended x, Extended y) {
return (!isunordered(x,y)) && (x < y);
}
inline bool islessequal(Extended x, Extended y) {
return (!isunordered(x,y)) && (x <= y);
}
inline bool islessgreater(Extended x, Extended y) {
return (!isunordered(x,y)) && ((x < y) || (x > y));
}
extern const Extended ExtendedPositiveInfinity;
extern const Extended ExtendedNegativeInfinity;
extern const Extended ExtendedNaN;
#endif
#endif
}
// eval method is 0 for x87, sse, soft, that is float_t has size float and double_t has size double
#define FLT_EVAL_METHOD 0
// Endianity checking
#if __FLOAT_WORD_ORDER == 1234
#define LOW_WORD_IDX 0
#define HIGH_WORD_IDX 1
#elif __FLOAT_WORD_ORDER == 4321
#define LOW_WORD_IDX 1
#define HIGH_WORD_IDX 0
#else
#error unknown byte order
#endif
#ifdef Extended
// little endian
#if __FLOAT_WORD_ORDER == 1234
/// Little endian is fine, always the same offsets whatever the memory model
template<int N> struct ExtendedConverter {
// Sign and exponent
static inline SizedUnsignedInteger<16>::Type* sexpPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<16>::Type*>(reinterpret_cast<char*>(e)+8);
}
// Mantissa
static inline SizedUnsignedInteger<32>::Type* mPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<32>::Type*>(e);
}
};
// Big endian
#elif __FLOAT_WORD_ORDER == 4321
/// Extended is softfloat
template<> struct ExtendedConverter<10> {
// Sign and exponent
static inline SizedUnsignedInteger<16>::Type* sexpPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<16>::Type*>(e);
}
// Mantissa
static inline SizedUnsignedInteger<32>::Type* mPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<32>::Type*>(reinterpret_cast<char*>(e)+2);
}
};
/// Default gcc model for x87 - 32 bits (or with -m96bit-long-double)
template<> struct ExtendedConverter<12> {
// Sign and exponent
static inline SizedUnsignedInteger<16>::Type* sexpPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<16>::Type*>(reinterpret_cast<char*>(e)+2);
}
// Mantissa
static inline SizedUnsignedInteger<32>::Type* mPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<32>::Type*>(reinterpret_cast<char*>(e)+4);
}
};
/// Default gcc model for x87 - 64 bits (or with -m128bit-long-double)
template<> struct ExtendedConverter<16> {
// Sign and exponent
static inline SizedUnsignedInteger<16>::Type* sexpPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<16>::Type*>(reinterpret_cast<char*>(e)+6);
}
// Mantissa
static inline SizedUnsignedInteger<32>::Type* mPtr(Extended* e) {
return reinterpret_cast<SizedUnsignedInteger<32>::Type*>(reinterpret_cast<char*>(e)+8);
}
};
#else
#error unknown byte order
#endif
#endif
// Simple
#ifdef LIBM_COMPILING_FLT32
// SSE is best case, always plain types
// X87 uses a wrapper for no denormals case
// First member of struct guaranted aligned at mem location => OK
// Idem for SoftFloat wrapper, though endianity has to be checked explicitly
#define SIMPLE_FROM_INT_PTR(x) *reinterpret_cast<Simple*>(x)
#define CONST_SIMPLE_FROM_INT_PTR(x) *reinterpret_cast<const Simple*>(x)
#define GET_FLOAT_WORD(i,d) \
do { \
Simple f = (d); \
(i) = *reinterpret_cast<streflop::uint32_t*>(&f); \
} while (0)
/* Set a float from a 32 bit int. */
#define SET_FLOAT_WORD(d,i) \
do { \
int ii = (i); \
(d) = *reinterpret_cast<Simple*>(&ii); \
} while (0)
#endif
// Double
#ifdef LIBM_COMPILING_DBL64
#define DOUBLE_FROM_INT_PTR(x) *reinterpret_cast<Double*>(x)
#define CONST_DOUBLE_FROM_INT_PTR(x) *reinterpret_cast<const Double*>(x)
#define EXTRACT_WORDS(ix0,ix1,d) \
do { \
Double f = (d); \
(ix0) = reinterpret_cast<streflop::uint32_t*>(&f)[HIGH_WORD_IDX]; \
(ix1) = reinterpret_cast<streflop::uint32_t*>(&f)[LOW_WORD_IDX]; \
} while (0)
#define GET_HIGH_WORD(i,d) \
do { \
Double f = (d); \
(i) = reinterpret_cast<streflop::uint32_t*>(&f)[HIGH_WORD_IDX]; \
} while (0)
#define GET_LOW_WORD(i,d) \
do { \
Double f = (d); \
(i) = reinterpret_cast<streflop::uint32_t*>(&f)[LOW_WORD_IDX]; \
} while (0)
#define INSERT_WORDS(d,ix0,ix1) \
do { \
Double f; \
reinterpret_cast<streflop::uint32_t*>(&f)[HIGH_WORD_IDX] = (ix0); \
reinterpret_cast<streflop::uint32_t*>(&f)[LOW_WORD_IDX] = (ix1); \
(d) = f; \
} while (0)
#define SET_HIGH_WORD(d,v) \
do { \
Double f = (d); \
reinterpret_cast<streflop::uint32_t*>(&f)[HIGH_WORD_IDX] = (v); \
(d) = f; \
} while (0)
#define SET_LOW_WORD(d,v) \
do { \
Double f = (d); \
reinterpret_cast<streflop::uint32_t*>(&f)[LOW_WORD_IDX] = (v); \
(d) = f; \
} while (0)
#endif
// Extended
#ifdef LIBM_COMPILING_LDBL96
#define EXTENDED_FROM_INT_PTR(x) *reinterpret_cast<Extended*>(x)
#define CONST_EXTENDED_FROM_INT_PTR(x) *reinterpret_cast<const Extended*>(x)
#define GET_LDOUBLE_WORDS(exp,ix0,ix1,d) \
do { \
Extended f = (d); \
(exp) = *ExtendedConverter<sizeof(Extended)>::sexpPtr(&f); \
(ix0) = ExtendedConverter<sizeof(Extended)>::mPtr(&f)[HIGH_WORD_IDX]; \
(ix1) = ExtendedConverter<sizeof(Extended)>::mPtr(&f)[LOW_WORD_IDX]; \
} while (0)
#define SET_LDOUBLE_WORDS(d,exp,ix0,ix1) \
do { \
Extended f; \
*ExtendedConverter<sizeof(Extended)>::sexpPtr(&f) = (exp); \
ExtendedConverter<sizeof(Extended)>::mPtr(&f)[HIGH_WORD_IDX] = (ix0); \
ExtendedConverter<sizeof(Extended)>::mPtr(&f)[LOW_WORD_IDX] = (ix1); \
(d) = f; \
} while (0)
#define GET_LDOUBLE_MSW(v,d) \
do { \
Extended f = (d); \
(v) = ExtendedConverter<sizeof(Extended)>::mPtr(&f)[HIGH_WORD_IDX]; \
} while (0)
#define SET_LDOUBLE_MSW(d,v) \
do { \
Extended f = (d); \
ExtendedConverter<sizeof(Extended)>::mPtr(&f)[HIGH_WORD_IDX] = (v); \
(d) = f; \
} while (0)
#define GET_LDOUBLE_EXP(exp,d) \
do { \
Extended f = (d); \
(exp) = *ExtendedConverter<sizeof(Extended)>::sexpPtr(&f); \
} while (0)
#define SET_LDOUBLE_EXP(d,exp) \
do { \
Extended f = (d); \
*ExtendedConverter<sizeof(Extended)>::sexpPtr(&f) = (exp); \
(d) = f; \
} while (0)
#endif
#ifdef _MSC_VER
// Even if MSVC != STDC, we still need to define it,
// otherwise MSVC chokes on the K&R style function headers.
# define __STDC__
#endif
#endif

99
source/shared_lib/include/streflop/streflop.h Normal file
View File

@@ -0,0 +1,99 @@
/*
streflop: STandalone REproducible FLOating-Point
Nicolas Brodu, 2006
Code released according to the GNU Lesser General Public License
Heavily relies on GNU Libm, itself depending on netlib fplibm, GNU MP, and IBM MP lib.
Uses SoftFloat too.
Please read the history and copyright information in the documentation provided with the source code
*/
#ifndef STREFLOP_H
#define STREFLOP_H
// protect against bad defines
#if defined(STREFLOP_SSE) && defined(STREFLOP_X87)
#error You have to define exactly one of STREFLOP_SSE STREFLOP_X87 STREFLOP_SOFT, but you defined both STREFLOP_SSE and STREFLOP_X87
#elif defined(STREFLOP_SSE) && defined(STREFLOP_SOFT)
#error You have to define exactly one of STREFLOP_SSE STREFLOP_X87 STREFLOP_SOFT, but you defined both STREFLOP_SSE and STREFLOP_SOFT
#elif defined(STREFLOP_X87) && defined(STREFLOP_SOFT)
#error You have to define exactly one of STREFLOP_SSE STREFLOP_X87 STREFLOP_SOFT, but you defined both STREFLOP_X87 and STREFLOP_SOFT
#elif !defined(STREFLOP_SSE) && !defined(STREFLOP_X87) && !defined(STREFLOP_SOFT)
#error You have to define exactly one of STREFLOP_SSE STREFLOP_X87 STREFLOP_SOFT, but you defined none
#endif
// First, define the numerical types
namespace streflop {
// Handle the 6 cells of the configuration array. See README.txt
#if defined(STREFLOP_SSE)
// SSE always uses native types, denormals are handled by FPU flags
typedef float Simple;
typedef double Double;
#undef Extended
#elif defined(STREFLOP_X87)
// X87 uses a wrapper for no denormals case
#if defined(STREFLOP_NO_DENORMALS)
#include "X87DenormalSquasher.h"
typedef X87DenormalSquasher<float> Simple;
typedef X87DenormalSquasher<double> Double;
typedef X87DenormalSquasher<long double> Extended;
#define Extended Extended
#else
// Use FPU flags for x87 with denormals
typedef float Simple;
typedef double Double;
typedef long double Extended;
#define Extended Extended
#endif
#elif defined(STREFLOP_SOFT) && !defined(STREFLOP_NO_DENORMALS)
// Use SoftFloat wrapper
#include "SoftFloatWrapper.h"
typedef SoftFloatWrapper<32> Simple;
typedef SoftFloatWrapper<64> Double;
typedef SoftFloatWrapper<96> Extended;
#define Extended Extended
#else
#error STREFLOP: Invalid combination or unknown FPU type.
#endif
}
#if defined(STREFLOP_SSE) && defined(_MSC_VER)
// MSVC will not compile without the long double variants defined, so we either have to declare Extended:
// typedef long double Extended;
// #define Extended Extended
// or write some wrapper macros:
#define atan2l(a, b) atan2((double)a, (double)b)
#define cosl(a) cos((double)a)
#define expl(a) exp((double)a)
#define ldexpl(a, b) ldexp((double)a, (double)b)
#define logl(a) log((double)a)
#define powl(a,b) pow((double)a, (double)b)
#define sinl(a) sin((double)a)
#define sqrtl(a) sqrt((double)a)
#define tanl(a) tan((double)a)
#define frexpl(a, b) frexp((double)a, b)
// The wrappers are possibly the better choice for sync reasons.
#endif
// Include the FPU settings file, so the user can initialize the library
#include "FPUSettings.h"
// Now that types are defined, include the Math.h file for the prototypes
#include "SMath.h"
// And now that math functions are defined, include the random numbers
#include "Random.h"
#endif

29
source/shared_lib/include/streflop/streflopC.h Normal file
View File

@@ -0,0 +1,29 @@
/* Copyright (C) 2009 Tobi Vollebregt */
/*
Serves as a C compatible interface to the most basic streflop functions.
*/
#ifndef STREFLOP_C_H
#define STREFLOP_C_H
#ifdef __cplusplus
extern "C" {
#endif
/// Initializes the FPU to single precision
void streflop_init_Simple();
/// Initializes the FPU to double precision
void streflop_init_Double();
#if defined(Extended)
/// Initializes the FPU to extended precision
void streflop_init_Extended();
#endif // defined(Extended)
#ifdef __cplusplus
} // extern "C"
#endif
#endif // STREFLOP_C_H

65
source/shared_lib/include/streflop/streflop_cond.h Normal file
View File

@@ -0,0 +1,65 @@
/* Copyright (C) 2008 Tobi Vollebregt */
/* Conditionally include streflop depending on STREFLOP_* #defines:
If one of those is present, #include "streflop.h", otherwise #include <math.h>
When faced with ambiguous call errors with e.g. fabs, use math::function.
Add it to math namespace if it doesn't exist there yet. */
#ifndef STREFLOP_COND_H
#define STREFLOP_COND_H
#if defined(STREFLOP_X87) || defined(STREFLOP_SSE) || defined(STREFLOP_SOFT)
#include "streflop.h"
using namespace streflop;
namespace math {
using namespace streflop;
}
#else
#include <cmath>
namespace math {
using std::fabs;
// using std::sqrt;
using std::sin;
using std::cos;
using std::sinh;
using std::cosh;
using std::tan;
using std::tanh;
using std::asin;
using std::acos;
using std::atan;
using std::atan2;
using std::ceil;
using std::floor;
using std::fmod;
using std::pow;
using std::log;
using std::log10;
using std::exp;
using std::frexp;
using std::ldexp;
// the following are C99 functions -> not supported by VS C
#if !defined(_MSC_VER) || _MSC_VER < 1500
using std::isnan;
using std::isinf;
using std::isfinite;
#elif __cplusplus
template<typename T> inline bool isnan(T value) {
return value != value;
}
// requires include <limits>
template<typename T> inline bool isinf(T value) {
return std::numeric_limits<T>::has_infinity && value == std::numeric_limits<T>::infinity();
}
// requires include <limits>
template<typename T> inline bool isfinite(T value) {
return !isinf<T>(value);
}
#endif
}
#endif
#endif // STREFLOP_COND_H