aboutsummaryrefslogtreecommitdiff
path: root/Client/ThirdParty
diff options
context:
space:
mode:
Diffstat (limited to 'Client/ThirdParty')
-rw-r--r--Client/ThirdParty/fix32/fix32.cpp0
-rw-r--r--Client/ThirdParty/fix32/fix32.h7
-rw-r--r--Client/ThirdParty/fix32/fix32.hpp161
-rw-r--r--Client/ThirdParty/math-sll/LICENSE23
-rw-r--r--Client/ThirdParty/math-sll/Makefile69
-rw-r--r--Client/ThirdParty/math-sll/README51
-rw-r--r--Client/ThirdParty/math-sll/math-sll.c957
-rw-r--r--Client/ThirdParty/math-sll/math-sll.h770
8 files changed, 2031 insertions, 7 deletions
diff --git a/Client/ThirdParty/fix32/fix32.cpp b/Client/ThirdParty/fix32/fix32.cpp
deleted file mode 100644
index e69de29..0000000
--- a/Client/ThirdParty/fix32/fix32.cpp
+++ /dev/null
diff --git a/Client/ThirdParty/fix32/fix32.h b/Client/ThirdParty/fix32/fix32.h
deleted file mode 100644
index fb3d65c..0000000
--- a/Client/ThirdParty/fix32/fix32.h
+++ /dev/null
@@ -1,7 +0,0 @@
-#pragma once
-
-#include <stdint.h>
-
-// -2,147,483,648 to 2,147,483,647
-typedef int64_t fix32_t;
-
diff --git a/Client/ThirdParty/fix32/fix32.hpp b/Client/ThirdParty/fix32/fix32.hpp
new file mode 100644
index 0000000..7999ad4
--- /dev/null
+++ b/Client/ThirdParty/fix32/fix32.hpp
@@ -0,0 +1,161 @@
+#pragma once
+
+#include <stdint.h>
+extern "C" {
+#include "math-sll/math-sll.h"
+}
+
+typedef sll fixed32_t;
+
+// Q32.32
+class Fix32
+{
+public:
+ sll value;
+
+ inline Fix32() { value = 0; }
+ inline Fix32(const Fix32 &inValue) { value = inValue.value; }
+ inline Fix32(const float inValue) { value = dbl2sll(inValue); }
+ inline Fix32(const double inValue) { value = dbl2sll(inValue); }
+ inline Fix32(const int32_t inValue) { value = int2sll(inValue); }
+ inline Fix32(const sll inValue) { value = inValue; }
+
+ inline operator sll() const { return value; }
+ inline operator double() const { return sll2dbl(value); }
+ inline operator float() const { return (float)sll2dbl(value); }
+ inline operator int32_t() const { return (int32_t)sll2int(value); }
+
+ inline Fix32 & operator=(const Fix32 &rhs) { value = rhs.value; return *this; }
+ inline Fix32 & operator=(const sll rhs) { value = rhs; return *this; }
+ inline Fix32 & operator=(const double rhs) { value = dbl2sll(rhs); return *this; }
+ inline Fix32 & operator=(const float rhs) { value = (float)dbl2sll(rhs); return *this; }
+ inline Fix32 & operator=(const int32_t rhs) { value = int2sll(rhs); return *this; }
+
+ inline Fix32 & operator+=(const Fix32 &rhs) { value += rhs.value; return *this; }
+ inline Fix32 & operator+=(const sll rhs) { value += rhs; return *this; }
+ inline Fix32 & operator+=(const double rhs) { value += dbl2sll(rhs); return *this; }
+ inline Fix32 & operator+=(const float rhs) { value += (float)dbl2sll(rhs); return *this; }
+ inline Fix32 & operator+=(const int32_t rhs) { value += int2sll(rhs); return *this; }
+
+ inline Fix32 & operator-=(const Fix32 &rhs) { value -= rhs.value; return *this; }
+ inline Fix32 & operator-=(const sll rhs) { value -= rhs; return *this; }
+ inline Fix32 & operator-=(const double rhs) { value -= dbl2sll(rhs); return *this; }
+ inline Fix32 & operator-=(const float rhs) { value -= (float)dbl2sll(rhs); return *this; }
+ inline Fix32 & operator-=(const int32_t rhs) { value -= int2sll(rhs); return *this; }
+
+ inline Fix32 & operator*=(const Fix32 &rhs) { value = sllmul(value, rhs.value); return *this; }
+ inline Fix32 & operator*=(const sll rhs) { value = sllmul(value, rhs); return *this; }
+ inline Fix32 & operator*=(const double rhs) { value = sllmul(value, dbl2sll(rhs)); return *this; }
+ inline Fix32 & operator*=(const float rhs) { value = sllmul(value, (float)dbl2sll(rhs)); return *this; }
+ inline Fix32 & operator*=(const int32_t rhs) { value *= rhs; return *this; }
+
+ inline Fix32 & operator/=(const Fix32 &rhs) { value = slldiv(value, rhs.value); return *this; }
+ inline Fix32 & operator/=(const sll rhs) { value = slldiv(value, rhs); return *this; }
+ inline Fix32 & operator/=(const double rhs) { value = slldiv(value, dbl2sll(rhs)); return *this; }
+ inline Fix32 & operator/=(const float rhs) { value = slldiv(value, (float)dbl2sll(rhs)); return *this; }
+ inline Fix32 & operator/=(const int32_t rhs) { value /= rhs; return *this; }
+
+ inline const Fix32 operator+(const Fix32 &other) const { Fix32 ret = *this; ret += other; return ret; }
+ inline const Fix32 operator+(const sll other) const { Fix32 ret = *this; ret += other; return ret; }
+ inline const Fix32 operator+(const double other) const { Fix32 ret = *this; ret += other; return ret; }
+ inline const Fix32 operator+(const float other) const { Fix32 ret = *this; ret += other; return ret; }
+ inline const Fix32 operator+(const int32_t other) const { Fix32 ret = *this; ret += other; return ret; }
+
+#ifndef FIXMATH_NO_OVERFLOW
+ inline const Fix32 sadd(const Fix32 &other) const { Fix32 ret = slladd(value, other.value); return ret; }
+ inline const Fix32 sadd(const sll other) const { Fix32 ret = slladd(value, other); return ret; }
+ inline const Fix32 sadd(const double other) const { Fix32 ret = slladd(value, dbl2sll(other)); return ret; }
+ inline const Fix32 sadd(const float other) const { Fix32 ret = slladd(value, (float)dbl2sll(other)); return ret; }
+ inline const Fix32 sadd(const int32_t other) const { Fix32 ret = slladd(value, int2sll(other)); return ret; }
+#endif
+
+ inline const Fix32 operator-(const Fix32 &other) const { Fix32 ret = *this; ret -= other; return ret; }
+ inline const Fix32 operator-(const sll other) const { Fix32 ret = *this; ret -= other; return ret; }
+ inline const Fix32 operator-(const double other) const { Fix32 ret = *this; ret -= other; return ret; }
+ inline const Fix32 operator-(const float other) const { Fix32 ret = *this; ret -= other; return ret; }
+ inline const Fix32 operator-(const int32_t other) const { Fix32 ret = *this; ret -= other; return ret; }
+
+ inline const Fix32 operator-() const { Fix32 ret = sllsub(0, value); return ret; }
+
+#ifndef FIXMATH_NO_OVERFLOW
+ inline const Fix32 ssub(const Fix32 &other) const { Fix32 ret = slladd(value, -other.value); return ret; }
+ inline const Fix32 ssub(const sll other) const { Fix32 ret = slladd(value, -other); return ret; }
+ inline const Fix32 ssub(const double other) const { Fix32 ret = slladd(value, -dbl2sll(other)); return ret; }
+ inline const Fix32 ssub(const float other) const { Fix32 ret = slladd(value, -(float)dbl2sll(other)); return ret; }
+ inline const Fix32 ssub(const int32_t other) const { Fix32 ret = slladd(value, -int2sll(other)); return ret; }
+#endif
+
+ inline const Fix32 operator*(const Fix32 &other) const { Fix32 ret = *this; ret *= other; return ret; }
+ inline const Fix32 operator*(const sll other) const { Fix32 ret = *this; ret *= other; return ret; }
+ inline const Fix32 operator*(const double other) const { Fix32 ret = *this; ret *= other; return ret; }
+ inline const Fix32 operator*(const float other) const { Fix32 ret = *this; ret *= other; return ret; }
+ inline const Fix32 operator*(const int32_t other) const { Fix32 ret = *this; ret *= other; return ret; }
+
+#ifndef FIXMATH_NO_OVERFLOW
+ inline const Fix32 smul(const Fix32 &other) const { Fix32 ret = sllmul(value, other.value); return ret; }
+ inline const Fix32 smul(const sll other) const { Fix32 ret = sllmul(value, other); return ret; }
+ inline const Fix32 smul(const double other) const { Fix32 ret = sllmul(value, dbl2sll(other)); return ret; }
+ inline const Fix32 smul(const float other) const { Fix32 ret = sllmul(value, (float)dbl2sll(other)); return ret; }
+ inline const Fix32 smul(const int32_t other) const { Fix32 ret = sllmul(value, int2sll(other)); return ret; }
+#endif
+
+ inline const Fix32 operator/(const Fix32 &other) const { Fix32 ret = *this; ret /= other; return ret; }
+ inline const Fix32 operator/(const sll other) const { Fix32 ret = *this; ret /= other; return ret; }
+ inline const Fix32 operator/(const double other) const { Fix32 ret = *this; ret /= other; return ret; }
+ inline const Fix32 operator/(const float other) const { Fix32 ret = *this; ret /= other; return ret; }
+ inline const Fix32 operator/(const int32_t other) const { Fix32 ret = *this; ret /= other; return ret; }
+
+#ifndef FIXMATH_NO_OVERFLOW
+ inline const Fix32 sdiv(const Fix32 &other) const { Fix32 ret = slldiv(value, other.value); return ret; }
+ inline const Fix32 sdiv(const sll other) const { Fix32 ret = slldiv(value, other); return ret; }
+ inline const Fix32 sdiv(const double other) const { Fix32 ret = slldiv(value, dbl2sll(other)); return ret; }
+ inline const Fix32 sdiv(const float other) const { Fix32 ret = slldiv(value, (float)dbl2sll(other)); return ret; }
+ inline const Fix32 sdiv(const int32_t other) const { Fix32 ret = slldiv(value, int2sll(other)); return ret; }
+#endif
+
+ inline int operator==(const Fix32 &other) const { return (value == other.value); }
+ inline int operator==(const sll other) const { return (value == other); }
+ inline int operator==(const double other) const { return (value == dbl2sll(other)); }
+ inline int operator==(const float other) const { return (value == (float)dbl2sll(other)); }
+ inline int operator==(const int32_t other) const { return (value == int2sll(other)); }
+
+ inline int operator!=(const Fix32 &other) const { return (value != other.value); }
+ inline int operator!=(const sll other) const { return (value != other); }
+ inline int operator!=(const double other) const { return (value != dbl2sll(other)); }
+ inline int operator!=(const float other) const { return (value != (float)dbl2sll(other)); }
+ inline int operator!=(const int32_t other) const { return (value != int2sll(other)); }
+
+ inline int operator<=(const Fix32 &other) const { return (value <= other.value); }
+ inline int operator<=(const sll other) const { return (value <= other); }
+ inline int operator<=(const double other) const { return (value <= dbl2sll(other)); }
+ inline int operator<=(const float other) const { return (value <= (float)dbl2sll(other)); }
+ inline int operator<=(const int32_t other) const { return (value <= int2sll(other)); }
+
+ inline int operator>=(const Fix32 &other) const { return (value >= other.value); }
+ inline int operator>=(const sll other) const { return (value >= other); }
+ inline int operator>=(const double other) const { return (value >= dbl2sll(other)); }
+ inline int operator>=(const float other) const { return (value >= (float)dbl2sll(other)); }
+ inline int operator>=(const int32_t other) const { return (value >= int2sll(other)); }
+
+ inline int operator< (const Fix32 &other) const { return (value < other.value); }
+ inline int operator< (const sll other) const { return (value < other); }
+ inline int operator< (const double other) const { return (value < dbl2sll(other)); }
+ inline int operator< (const float other) const { return (value < (float)dbl2sll(other)); }
+ inline int operator< (const int32_t other) const { return (value < int2sll(other)); }
+
+ inline int operator> (const Fix32 &other) const { return (value > other.value); }
+ inline int operator> (const sll other) const { return (value > other); }
+ inline int operator> (const double other) const { return (value > dbl2sll(other)); }
+ inline int operator> (const float other) const { return (value > (float)dbl2sll(other)); }
+ inline int operator> (const int32_t other) const { return (value > int2sll(other)); }
+
+ inline Fix32 sin() const { return Fix32(sllsin(value)); }
+ inline Fix32 cos() const { return Fix32(sllcos(value)); }
+ inline Fix32 tan() const { return Fix32(slltan(value)); }
+ inline Fix32 asin() const { return Fix32(sllasin(value)); }
+ inline Fix32 acos() const { return Fix32(sllacos(value)); }
+ inline Fix32 atan() const { return Fix32(sllatan(value)); }
+ //inline Fix32 atan2(const Fix32 &inY) const { return Fix32(fix16_atan2(value, inY.value)); }
+ inline Fix32 sqrt() const { return Fix32(sllsqrt(value)); }
+
+}; \ No newline at end of file
diff --git a/Client/ThirdParty/math-sll/LICENSE b/Client/ThirdParty/math-sll/LICENSE
new file mode 100644
index 0000000..4e683cf
--- /dev/null
+++ b/Client/ThirdParty/math-sll/LICENSE
@@ -0,0 +1,23 @@
+
+Licensed under the terms of the MIT license:
+
+Copyright (c) 2000,2002,2006,2012,2016 Andrew E. Mileski <andrewm@isoar.ca>
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to
+deal in the Software without restriction, including without limitation the
+rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
+sell copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The copyright notice, and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+DEALINGS IN THE SOFTWARE.
+
diff --git a/Client/ThirdParty/math-sll/Makefile b/Client/ThirdParty/math-sll/Makefile
new file mode 100644
index 0000000..b7b512f
--- /dev/null
+++ b/Client/ThirdParty/math-sll/Makefile
@@ -0,0 +1,69 @@
+#
+# Licensed under the terms of the MIT license:
+#
+# Copyright (c) 2000,2002,2006,2012,2016 Andrew E. Mileski <andrewm@isoar.ca>
+#
+# Permission is hereby granted, free of charge, to any person obtaining a copy
+# of this software and associated documentation files (the "Software"), to
+# deal in the Software without restriction, including without limitation the
+# rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
+# sell copies of the Software, and to permit persons to whom the Software is
+# furnished to do so, subject to the following conditions:
+#
+# The copyright notice, and this permission notice shall be included in all
+# copies or substantial portions of the Software.
+#
+# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+# FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+# DEALINGS IN THE SOFTWARE.
+#
+
+#
+# Installation directories
+#
+
+PREFIX := /usr/local
+INCDIR := $(PREFIX)/include
+LIBDIR := $(PREFIX)/lib
+
+#
+# Executables
+#
+
+AR := ar
+CC := gcc
+CFLAGS := -O2 -W -Wall
+INSTALL := install
+RANLIB := ranlib
+RM := rm -f
+STRIP := strip --strip-unneeded
+
+#
+# Recipes
+#
+
+OBJS := math-sll.o
+LIBS := math-sll.a
+
+.PHONY: all clean install
+
+all: $(LIBS)
+
+clean:
+ $(RM) $(LIBS) $(OBJS)
+
+install: $(LIBS) math-sll.h
+ $(INSTALL) -m a=rx,u+w math-sll.a $(LIBDIR)
+ $(INSTALL) -m a=r,u+w math-sll.h $(INCDIR)
+
+math-sll.o: math-sll.c math-sll.h
+
+math-sll.a: math-sll.o
+ $(STRIP) $<
+ $(AR) rcs $@ $<
+ $(RANLIB) $@
+
diff --git a/Client/ThirdParty/math-sll/README b/Client/ThirdParty/math-sll/README
new file mode 100644
index 0000000..e53ca36
--- /dev/null
+++ b/Client/ThirdParty/math-sll/README
@@ -0,0 +1,51 @@
+
+math-sll
+
+ A fixed point (32.32 bit) math library.
+
+ See math-sll.h for details.
+
+License
+
+ Licensed under the terms of the MIT license as of revision v1.19
+
+ See the LICENSE file for details.
+
+ Earlier revisions were licensed under the terms of GPL or LGPL.
+
+Installation
+
+ The source can be added as-is into your project, or a library can optionally
+ be built and installed, then statically-linked into your project.
+
+ To build:
+
+ make clean
+ make all
+
+ Optionally:
+
+ make install
+
+ The default installation path prefix is /usr/local
+
+ See the Makefile for details.
+
+Repository
+
+ A math-sll project is being maintained at:
+
+ http://github.com/aemileski/math-sll
+
+ You can also obtain a copy of the latest source code directly via GIT:
+
+ git clone git://github.com/aemileski/math-sll.git
+
+Contact
+
+ Feel free to contact me with questions, bug reports, patches, feature
+ requests, etc., or just to tell me about how you are using math-sll in
+ your project!
+
+ Andrew E. Mileski <andrewm@isoar.ca>
+
diff --git a/Client/ThirdParty/math-sll/math-sll.c b/Client/ThirdParty/math-sll/math-sll.c
new file mode 100644
index 0000000..ecc4238
--- /dev/null
+++ b/Client/ThirdParty/math-sll/math-sll.c
@@ -0,0 +1,957 @@
+/*
+ * Revision v1.24
+ *
+ * Credits
+ *
+ * Maintained, conceived, written, and fiddled with by:
+ *
+ * Andrew E. Mileski <andrewm@isoar.ca>
+ *
+ * Other source code contributors:
+ *
+ * Kevin Rockel
+ * Kevin Michael Woley
+ * Mark Anthony Lisher
+ * Nicolas Pitre
+ * Anonymous
+ *
+ * License
+ *
+ * Licensed under the terms of the MIT license:
+ *
+ * Copyright (c) 2000,2002,2006,2012,2016 Andrew E. Mileski <andrewm@isoar.ca>
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to
+ * deal in the Software without restriction, including without limitation the
+ * rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
+ * sell copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The copyright notice, and this permission notice shall be included in all
+ * copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+ * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+ * DEALINGS IN THE SOFTWARE.
+ */
+
+/* See header for full details */
+#include "math-sll.h"
+
+/*
+ * Local prototypes
+ */
+
+static sll _sllcos(sll x);
+static sll _sllsin(sll x);
+
+static sll _sllexp(sll x);
+
+/*
+ * Unpack IEEE 754 floating point double format into fixed point sll format
+ *
+ * Description
+ *
+ * IEEE 754 specifies the binary64 type ("double" in C) as having:
+ *
+ * 1 bit sign
+ * 11 bit exponent
+ * 53 bit significand
+ *
+ * The first bit of the significand is an implied 1 which is not stored.
+ * The decimal would be to the right of that implied 1, or to the left of
+ * the stored significand.
+ *
+ * The exponent is unsigned, and biased with an offset of 1023.
+ *
+ * The IEEE 754 standard does not specify endianess, but the endian used is
+ * traditionally the same endian that the processor uses.
+ */
+
+sll dbl2sll(double dbl)
+{
+ union {
+ double d;
+ unsigned u[2];
+ ull ull;
+ sll sll;
+ } in, retval;
+ register unsigned exp;
+
+ /* Move into memory as args might be passed in regs */
+ in.d = dbl;
+
+#if defined(BROKEN_IEEE754_DOUBLE)
+
+ exp = in.u[0];
+ in.u[0] = in.u[1];
+ in.u[1] = exp;
+
+#endif /* defined(BROKEN_IEEE754_DOUBLE) */
+
+ /* Leading 1 is assumed by IEEE */
+ retval.u[1] = 0x40000000;
+
+ /* Unpack the mantissa into the unsigned long */
+ retval.u[1] |= (in.u[1] << 10) & 0x3ffffc00;
+ retval.u[1] |= (in.u[0] >> 22) & 0x000003ff;
+ retval.u[0] = in.u[0] << 10;
+
+ /* Extract the exponent and align the decimals */
+ exp = (in.u[1] >> 20) & 0x7ff;
+ if (exp)
+ /* IEEE 754 decimal begins at right of bit position 30 */
+ retval.ull >>= (1023 + 30) - exp;
+ else
+ return 0L;
+
+ /* Negate if negative flag set */
+ if (in.u[1] & 0x80000000)
+ retval.sll = _sllneg(retval.sll);
+
+ return retval.sll;
+}
+
+/*
+ * Pack fixed point sll format into IEEE 754 floating point double format
+ *
+ * Description
+ *
+ * IEEE 754 specifies the binary64 type ("double" in C) as having:
+ *
+ * 1 bit sign
+ * 11 bit exponent
+ * 53 bit significand
+ *
+ * The first bit of the significand is an implied 1 which is not stored.
+ * The decimal would be to the right of that implied 1, or to the left of
+ * the stored significand.
+ *
+ * The exponent is unsigned, and biased with an offset of 1023.
+ *
+ * The IEEE 754 standard does not specify endianess, but the endian used is
+ * traditionally the same endian that the processor uses.
+ */
+
+double sll2dbl(sll s)
+{
+ union {
+ double d;
+ unsigned u[2];
+ ull ull;
+ sll sll;
+ } in, retval;
+ register unsigned exp;
+ register unsigned flag;
+
+ if (s == 0)
+ return 0.0;
+
+ /* Move into memory as args might be passed in regs */
+ in.sll = s;
+
+ /* Handle the negative flag */
+ if (in.sll < 1) {
+ flag = 0x80000000;
+ in.ull = _sllneg(in.sll);
+ } else
+ flag = 0x00000000;
+
+ /*
+ * Normalize
+ *
+ * IEEE 754 decimal-point begins at right of bit position 30
+ */
+ for (exp = (1023 + 30); in.ull && (in.u[1] & 0x80000000) == 0; exp--) {
+ in.ull <<= 1;
+ }
+ in.ull <<= 1;
+ exp++;
+ in.ull >>= 12;
+ retval.ull = in.ull;
+ retval.u[1] |= flag | (exp << 20);
+
+#if defined(BROKEN_IEEE754_DOUBLE)
+
+ exp = retval.u[0];
+ retval.u[0] = retval.u[1];
+ retval.u[1] = exp;
+
+#endif /* defined(BROKEN_IEEE754_DOUBLE) */
+
+ return retval.d;
+}
+
+/*
+ * Multiply two sll values
+ *
+ * Description
+ *
+ * When multiplying two 64 bit sll numbers, the result is 128 bits, but there
+ * is only room for a 64 bit result with sll!
+ *
+ * The 128 bit result has 64 bits on either side of the decimal, so 32 bits
+ * of overflow to the left of the decimal, and 32 bits of underflow to the
+ * right of the decmial.
+ *
+ * 32.32 * 32.32 = 64.64 = overflow(32) + 32.32 + underflow(32)
+ *
+ * However, a "long long" multiply has 64 bits of overflow to the left of the
+ * decimal, resulting in the entire integer part being lost!
+ *
+ * 32.32 * 32.32 = 64.64 = overflow(64) + .64
+ *
+ * Hence a custom multiply routine is required, to preserve the parts
+ * of the result that sll needs.
+ *
+ * Consider two sll numbers, x and y:
+ *
+ * Let x = x_hi * 2^0 + x_lo * 2^(-32)
+ * Let y = y_hi * 2^0 + y_lo * 2^(-32)
+ *
+ * Where:
+ *
+ * *_hi is the signed 32 bit integer part to the left of the decimal
+ * *_lo is the unsigned 32 bit fractional part to the right of the decimal
+ *
+ * x * y = (x_hi * 2^0 + x_lo * 2^(-32))
+ * * (y_hi * 2^0 + y_lo * 2^(-32))
+ *
+ * Expanding the terms, we get:
+ *
+ * = x_hi * y_hi * 2^0 + x_hi * y_lo * 2^(-32)
+ * + x_lo * y_hi * 2^(-32) + x_lo * y_lo * 2^(-64)
+ *
+ * Grouping by powers of 2, we get:
+ *
+ * (x_hi * y_hi) * 2^0
+ * We only need the low 32 bits of this term, as the rest is overflow
+ *
+ * (x_hi * y_lo + x_lo * y_hi) * 2^-32
+ * We need all bits of this term
+ *
+ * x_lo * y_lo * 2^-64
+ * We only need the high 32 bits of this term, as the rest is underflow
+ */
+
+sll sllmul(sll x, sll y)
+{
+ register unsigned int x_lo;
+ register signed int x_hi;
+
+ register unsigned int y_lo;
+ register signed int y_hi;
+
+ x_hi = (signed int) ((ull) x >> 32); // Discard lower 32 bits
+ x_lo = (unsigned int) x; // Discard upper 32 bits
+
+ y_hi = (signed int) ((ull) y >> 32); // Discard lower 32 bits
+ y_lo = (unsigned int) y; // Discard upper 32 bits
+
+ return (sll) (
+ ((ull) (x_hi * y_hi) << 32)
+ + ((ull) x_hi * y_lo + x_lo * (ull) y_hi)
+ + (((ull) x_lo * y_lo) >> 32)
+ );
+}
+
+/*
+ * Calculate cos x where -pi/4 <= x <= pi/4
+ *
+ * Description
+ *
+ * cos x = 1 - x^2 / 2! + x^4 / 4! - ... + x^(2N) / (2N)!
+ * Note that (pi/4)^12 / 12! < 2^-32 which is the smallest possible number.
+ *
+ * cos x = t0 + t1 + t2 + t3 + t4 + t5 + t6
+ *
+ * Consider only the factorials:
+ * f0 = 0! = 1
+ * f1 = 2! = 2 * 1 * f0 = 2 * f0
+ * f2 = 4! = 4 * 3 * f1 = 12 * f1
+ * f3 = 6! = 6 * 5 * f2 = 30 * f2
+ * f4 = 8! = 8 * 7 * f3 = 56 * f3
+ * f5 = 10! = 10 * 9 * f4 = 90 * f4
+ * f6 = 12! = 12 * 11 * f5 = 132 * f5
+ *
+ * Now consider each term of the series:
+ * t0 = 1
+ * t1 = -t0 * x^2 / f1 = -t0 * x^2 * CONST_1_2
+ * t2 = -t1 * x^2 / f2 = -t1 * x^2 * CONST_1_12
+ * t3 = -t2 * x^2 / f3 = -t2 * x^2 * CONST_1_30
+ * t4 = -t3 * x^2 / f4 = -t3 * x^2 * CONST_1_56
+ * t5 = -t4 * x^2 / f5 = -t4 * x^2 * CONST_1_90
+ * t6 = -t5 * x^2 / f6 = -t5 * x^2 * CONST_1_132
+ */
+
+sll _sllcos(sll x)
+{
+ sll retval;
+ sll x2;
+
+ x2 = sllmul(x, x);
+
+ retval = _sllsub(CONST_1, sllmul(x2, CONST_1_132));
+ retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_90));
+ retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_56));
+ retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_30));
+ retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_12));
+ retval = _sllsub(CONST_1, slldiv2(sllmul(x2, retval)));
+
+ return retval;
+}
+
+/*
+ * Calculate sin x where -pi/4 <= x <= pi/4
+ *
+ * Description
+ *
+ * sin x = x - x^3 / 3! + x^5 / 5! - ... + x^(2N+1) / (2N+1)!
+ * Note that (pi/4)^13 / 13! < 2^-32 which is the smallest possible number.
+ *
+ * sin x = t0 + t1 + t2 + t3 + t4 + t5 + t6
+ *
+ * Consider only the factorials:
+ * f0 = 0! = 1
+ * f1 = 3! = 3 * 2 * f0 = 6 * f0
+ * f2 = 5! = 5 * 4 * f1 = 20 * f1
+ * f3 = 7! = 7 * 6 * f2 = 42 * f2
+ * f4 = 9! = 9 * 8 * f3 = 72 * f3
+ * f5 = 11! = 11 * 10 * f4 = 110 * f4
+ * f6 = 13! = 13 * 12 * f5 = 156 * f5
+ *
+ * Now consider each term of the series:
+ * t0 = 1
+ * t1 = -t0 * x^2 / 6 = -t0 * x^2 * CONST_1_6
+ * t2 = -t1 * x^2 / 20 = -t1 * x^2 * CONST_1_20
+ * t3 = -t2 * x^2 / 42 = -t2 * x^2 * CONST_1_42
+ * t4 = -t3 * x^2 / 72 = -t3 * x^2 * CONST_1_72
+ * t5 = -t4 * x^2 / 110 = -t4 * x^2 * CONST_1_110
+ * t6 = -t5 * x^2 / 156 = -t5 * x^2 * CONST_1_156
+ */
+
+sll _sllsin(sll x)
+{
+ sll retval;
+ sll x2;
+
+ x2 = sllmul(x, x);
+
+ retval = _sllsub(x, sllmul(x2, CONST_1_156));
+ retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_110));
+ retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_72));
+ retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_42));
+ retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_20));
+ retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_6));
+
+ return retval;
+}
+
+/*
+ * Calculate cos x for any value of x, by quadrant
+ */
+
+sll sllcos(sll x)
+{
+ int i;
+ sll retval;
+
+ /* Calculate cos (x - i * pi/2), where -pi/4 <= x - i * pi/2 <= pi/4 */
+ i = _sll2int(_slladd(sllmul(x, CONST_2_PI), CONST_1_2));
+ x = _sllsub(x, sllmul(_int2sll(i), CONST_PI_2));
+
+ /* Locate the quadrant */
+ switch (i & 3) {
+ default:
+ case 0:
+ retval = _sllcos(x);
+ break;
+ case 1:
+ retval = sllneg(_sllsin(x));
+ break;
+ case 2:
+ retval = sllneg(_sllcos(x));
+ break;
+ case 3:
+ retval = _sllsin(x);
+ break;
+ }
+
+ return retval;
+}
+
+/*
+ * Calculate sin x for any value of x, by quadrant
+ */
+
+sll sllsin(sll x)
+{
+ int i;
+ sll retval;
+
+ /* Calculate sin (x - n * pi/2), where -pi/4 <= x - i * pi/2 <= pi/4 */
+ i = _sll2int(_slladd(sllmul(x, CONST_2_PI), CONST_1_2));
+ x = _sllsub(x, sllmul(_int2sll(i), CONST_PI_2));
+
+ /* Locate the quadrant */
+ switch (i & 3) {
+ default:
+ case 0:
+ retval = _sllsin(x);
+ break;
+ case 1:
+ retval = _sllcos(x);
+ break;
+ case 2:
+ retval = sllneg(_sllsin(x));
+ break;
+ case 3:
+ retval = sllneg(_sllcos(x));
+ break;
+ }
+
+ return retval;
+}
+
+/*
+ * Calculate tan x for any value of x, by quadrant
+ */
+
+sll slltan(sll x)
+{
+ int i;
+ sll retval;
+
+ i = _sll2int(_slladd(sllmul(x, CONST_2_PI), CONST_1_2));
+ x = _sllsub(x, sllmul(_int2sll(i), CONST_PI_2));
+
+ /* Locate the quadrant */
+ switch (i & 3) {
+ default:
+ case 0:
+ case 2:
+ retval = slldiv(_sllsin(x), _sllcos(x));
+ break;
+ case 1:
+ case 3:
+ retval = _sllneg(slldiv(_sllcos(x), _sllsin(x)));
+ break;
+ }
+
+ return retval;
+}
+
+/*
+ *
+ * Calculate asin x, where |x| <= 1
+ *
+ * Description
+ *
+ * asin x = SUM[n=0,) C(2 * n, n) * x^(2 * n + 1) / (4^n * (2 * n + 1)), |x| <= 1
+ *
+ * where C(n, r) = nCr = n! / (r! * (n - r)!)
+ *
+ * Using a two term approximation:
+ * [1] a = x + x^3 / 6
+ *
+ * Results in:
+ * asin x = a + D
+ * where D is the difference from the exact result
+ *
+ * Letting D = asin d results in:
+ * [2] asin x = a + asin d
+ *
+ * Re-arranging:
+ * asin x - a = asin d
+ *
+ * Applying sin to both sides:
+ * sin (asin x - a) = sin asin d
+ * sin (asin x - a) = d
+ * d = sin (asin x - a)
+ *
+ * Applying the standard identity:
+ * sin (u - v) = sin u * cos v - cos u * sin v
+ *
+ * Results in:
+ * d = sin asin x * cos a - cos asin x * sin a
+ * d = x * cos a - cos asin x * sin a
+ *
+ * Applying the standard identity:
+ * cos asin u = (1 - u^2)^(1 / 2)
+ *
+ * Results in:
+ * [3] d = x * cos a - (1 - x^2)^(1 / 2) * sin a
+ *
+ * Putting the pieces together:
+ * [1] a = x + x^3 / 6
+ * [3] d = x * cos a - (1 - x^2)^(1 / 2) * sin a
+ * [2] asin x = a + asin d
+ *
+ * The worst case is x = 1.0 which converges after 2 iterations.
+ */
+
+sll sllasin(sll x)
+{
+ int left_side;
+ sll a;
+ sll retval;
+
+ /* asin -x = -asin x */
+ if ((left_side = x < 0))
+ x = _sllneg(x);
+
+ /* Out-of-range */
+ if (x > CONST_1)
+ return 0;
+
+ /* Initial approximate */
+ a = sllmul(x, _slladd(CONST_1, sllmul(x, sllmul(x, CONST_1_6))));
+ retval = a;
+
+ /* First iteration */
+ x = _sllsub(sllmul(x, sllcos(a)), sllmul(sllsqrt(_sllsub(CONST_1, sllmul(x, x))), sllsin(a)));
+ a = sllmul(x, _slladd(CONST_1, sllmul(x, sllmul(x, CONST_1_6))));
+ retval = _slladd(retval, a);
+
+ /* Second iteration */
+ x = _sllsub(sllmul(x, sllcos(a)), sllmul(sllsqrt(_sllsub(CONST_1, sllmul(x, x))), sllsin(a)));
+ a = sllmul(x, _slladd(CONST_1, sllmul(x, sllmul(x, CONST_1_6))));
+ retval = _slladd(retval, a);
+
+ /* Negate result if necessary */
+ return (left_side ? _sllneg(retval): retval);
+}
+
+/*
+ * Calculate atan x
+ *
+ * Description
+ *
+ * atan x = SUM[n=0,) (-1)^n * x^(2 * n + 1) / (2 * n + 1), |x| <= 1
+ *
+ * Using a two term approximation:
+ * [1] a = x - x^3 / 3
+ *
+ * Results in:
+ * atan x = a + D
+ * where D is the difference from the exact result
+ *
+ * Letting D = atan d results in:
+ * [2] atan x = a + atan d
+ *
+ * Re-arranging:
+ * atan x - a = atan d
+ *
+ * Applying tan to both sides:
+ * tan (atan x - a) = tan atan d
+ * tan (atan x - a) = d
+ * d = tan (atan x - a)
+ *
+ * Applying the standard identity:
+ * tan (u - v) = (tan u - tan v) / (1 + tan u * tan v)
+ *
+ * Results in:
+ * d = tan (atan x - a) = (tan atan x - tan a) / (1 + tan atan x * tan a)
+ * d = tan (atan x - a) = (x - tan a) / (1 + x * tan a)
+ *
+ * Let:
+ * [3] t = tan a
+ *
+ * Results in:
+ * [4] d = (x - t) / (1 + x * t)
+ *
+ * So putting the pieces together:
+ * [1] a = x - x^3 / 3
+ * [3] t = tan a
+ * [4] d = (x - t) / (1 + x * t)
+ * [2] atan x = a + atan d
+ * atan x = a + atan ((x - t) / (1 + x * t))
+ *
+ * The worst case is x = 1.0 which converges after 2 iterations.
+ */
+
+sll sllatan(sll x)
+{
+ int side;
+ sll a;
+ sll t;
+ sll retval;
+
+
+ if (x < CONST_1) {
+
+ /* Left: if (x < -1) then atan x = pi / 2 + atan 1 / x */
+ side = -1;
+ x = sllinv(x);
+
+ } else if (x > CONST_1) {
+
+ /* Right: if (x > 1) then atan x = pi / 2 - atan 1 / x */
+ side = 1;
+ x = sllinv(x);
+
+ } else {
+ /* Middle: -1 <= x <= 1 */
+ side = 0;
+ }
+
+ /* Initial approximate */
+ a = sllmul(x, _sllsub(CONST_1, sllmul(x, sllmul(x, CONST_1_3))));
+ retval = a;
+
+ /* First iteration */
+ t = _slldiv(_sllsin(a), _sllcos(a));
+ x = _slldiv(_sllsub(x, t), _slladd(CONST_1, sllmul(x, t)));
+ a = sllmul(x, _sllsub(CONST_1, sllmul(x, sllmul(x, CONST_1_3))));
+ retval = _slladd(retval, a);
+
+ /* Second iteration */
+ t = _slldiv(_sllsin(a), _sllcos(a));
+ x = _slldiv(_sllsub(x, t), _slladd(CONST_1, sllmul(x, t)));
+ a = sllmul(x, _sllsub(CONST_1, sllmul(x, sllmul(x, CONST_1_3))));
+ retval = _slladd(retval, a);
+
+ if (side == -1) {
+
+ /* Left: if (x < -1) then atan x = pi / 2 + atan 1 / x */
+ retval = _slladd(CONST_PI_2, retval);
+
+ } else if (side == 1) {
+
+ /* Right: if (x > 1) then atan x = pi / 2 - atan 1 / x */
+ retval = _sllsub(CONST_PI_2, retval);
+ }
+
+ return retval;
+}
+
+/*
+ * Calculate e^x where -0.5 <= x <= 0.5
+ *
+ * Description:
+ * e^x = x^0 / 0! + x^1 / 1! + ... + x^N / N!
+ * Note that 0.5^11 / 11! < 2^-32 which is the smallest possible number.
+ */
+
+sll _sllexp(sll x)
+{
+ sll retval;
+
+ retval = CONST_1;
+
+ retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_11)));
+ retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_10)));
+ retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_9)));
+ retval = _slladd(CONST_1, sllmul(retval, slldiv2n(x, 3)));
+ retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_7)));
+ retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_6)));
+ retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_5)));
+ retval = _slladd(CONST_1, sllmul(retval, slldiv4(x)));
+ retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_3)));
+ retval = _slladd(CONST_1, sllmul(retval, slldiv2(x)));
+ retval = _slladd(CONST_1, sllmul(retval, x));
+
+ return retval;
+}
+
+/*
+ * Calculate e^x for any value of x
+ */
+
+sll sllexp(sll x)
+{
+ int i;
+ sll e;
+ sll retval;
+
+ e = CONST_E;
+
+ /* -0.5 <= x <= 0.5 */
+ i = _sll2int(_slladd(x, CONST_1_2));
+ retval = _sllexp(_sllsub(x, _int2sll(i)));
+
+ /* i >= 0 */
+ if (i < 0) {
+ i = -i;
+ e = CONST_1_E;
+ }
+
+ /* Scale the result */
+ for (; i; i >>= 1) {
+ if (i & 1)
+ retval = sllmul(retval, e);
+ e = sllmul(e, e);
+ }
+
+ return retval;
+}
+
+/*
+ * Calculate natural logarithm using Netwton-Raphson method
+ */
+
+sll slllog(sll x)
+{
+ sll x1;
+ sll ln;
+
+ ln = 0;
+
+ /* Scale: e^(-1/2) <= x <= e^(1/2) */
+ while (x < CONST_1_SQRTE) {
+ ln = _sllsub(ln, CONST_1);
+ x = sllmul(x, CONST_E);
+ }
+ while (x > CONST_SQRTE) {
+ ln = _slladd(ln, CONST_1);
+ x = sllmul(x, CONST_1_E);
+ }
+
+ /* First iteration */
+ x1 = sllmul(_sllsub(x, CONST_1), slldiv2(_sllsub(x, CONST_3)));
+ ln = _sllsub(ln, x1);
+ x = sllmul(x, _sllexp(x1));
+
+ /* Second iteration */
+ x1 = sllmul(_sllsub(x, CONST_1), slldiv2(_sllsub(x, CONST_3)));
+ ln = _sllsub(ln, x1);
+ x = sllmul(x, _sllexp(x1));
+
+ /* Third iteration */
+ x1 = sllmul(_sllsub(x, CONST_1), slldiv2(_sllsub(x, CONST_3)));
+ ln = _sllsub(ln, x1);
+
+ return ln;
+}
+
+/*
+ * Calculate the inverse for non-zero values
+ */
+
+sll sllinv(sll x)
+{
+ int sgn;
+ sll u;
+ ull s;
+
+ /* Use positive numbers, or the approximation won't work */
+ if (x < CONST_0) {
+ x = _sllneg(x);
+ sgn = 1;
+ } else {
+ sgn = 0;
+ }
+
+ /* Starting-point (gets shifted right to become positive) */
+ s = -1;
+
+ /* An approximation - must be larger than the actual value */
+ for (u = x; u; u = ((ull) u) >> 1)
+ s >>= 1;
+
+ /* Newton's Method */
+ u = sllmul(s, _sllsub(CONST_2, sllmul(x, s)));
+ u = sllmul(u, _sllsub(CONST_2, sllmul(x, u)));
+ u = sllmul(u, _sllsub(CONST_2, sllmul(x, u)));
+ u = sllmul(u, _sllsub(CONST_2, sllmul(x, u)));
+ u = sllmul(u, _sllsub(CONST_2, sllmul(x, u)));
+ u = sllmul(u, _sllsub(CONST_2, sllmul(x, u)));
+
+ return ((sgn) ? _sllneg(u): u);
+}
+
+/*
+ * Calculate x^y
+ *
+ * Description
+ *
+ * The standard identity:
+ * ln x^y = y * log x
+ *
+ * Raising e to the power of either sides:
+ * e^(ln x^y) = e^(y * log x)
+ *
+ * Which simplifies to:
+ * x^y = e^(y * ln x)
+ */
+
+sll sllpow(sll x, sll y)
+{
+ if (y == CONST_0)
+ return CONST_1;
+ if (y == CONST_1)
+ return x;
+ if (y == CONST_2)
+ return sllmul(x, x);
+
+ return sllexp(sllmul(y, slllog(x)));
+}
+
+/*
+ * Calculate the square-root
+ *
+ * Description
+ *
+ * Consider a parabola centered on the y-axis:
+ * y = a * x^2 + b
+ *
+ * Has zeros (y = 0) located at:
+ * a * x^2 + b = 0
+ * a * x^2 = -b
+ * x^2 = -b / a
+ * x = +- (-b / a)^(1 / 2)
+ *
+ * Letting a = 1 and b = -X results in:
+ * y = x^2 - X
+ * x = +- X^(1 / 2)
+ *
+ * Which is a convenient result, since we want to find the square root of X,
+ * and we can * use Newton's Method to find the zeros of any f(x):
+ * xn = x - f(x) / f'(x)
+ *
+ * Applying Newton's Method to our parabola:
+ * f(x) = x^2 - X
+ * xn = x - (x^2 - X) / (2 * x)
+ * xn = x - (x - X / x) / 2
+ *
+ * To make this converge quickly, we scale X so that:
+ * X = 4^N * z
+ *
+ * Taking the roots of both sides
+ * X^(1 / 2) = (4^n * z)^(1 / 2)
+ * X^(1 / 2) = 2^n * z^(1 / 2)
+ *
+ * Letting N = 2^n results in:
+ * x^(1 / 2) = N * z^(1 / 2)
+ *
+ * We want this to converge to the positive root, so we must start at a point
+ * 0 < start <= x^(1 / 2)
+ * or
+ * x^(1/2) <= start <= infinity
+ *
+ * Since:
+ * (1/2)^(1/2) = 0.707
+ * 2^(1/2) = 1.414
+ *
+ * A good choice is 1 which lies in the middle, and takes 4 iterations to
+ * converge from either extreme.
+ */
+
+sll sllsqrt(sll x)
+{
+ sll n;
+ sll xn;
+
+ /* Quick solutions for the simple cases */
+ if (x <= CONST_0 || x == CONST_1)
+ return x;
+
+ /* Start with a scaling factor of 1 */
+ n = CONST_1;
+
+ /* Scale x so that 0.5 <= x < 2 */
+ while (x >= CONST_2) {
+ x = slldiv4(x);
+ n = sllmul2(n);
+ }
+ while (x < CONST_1_2) {
+ x = sllmul4(x);
+ n = slldiv2(n);
+ }
+
+ /* Simple solution if x = 4^n */
+ if (x == CONST_1)
+ return n;
+
+ /* The starting point */
+ xn = CONST_1;
+
+ /* Four iterations will be enough */
+ xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));
+ xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));
+ xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));
+ xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));
+
+ /* Scale the result */
+ return sllmul(n, xn);
+}
+
+sll slld2dsqrt(sll num)
+{
+ if (num <= CONST_0 || num == CONST_1)
+ {
+ return num;
+ }
+ sll result = 0;
+ sll bit;
+
+ // Many numbers will be less than 15, so
+ // this gives a good balance between time spent
+ // in if vs. time spent in the while loop
+ // when searching for the starting value.
+ if (num & 0xFFFFFF0000000000LL)
+ bit = 0x4000000000000000LL; //(sll)1 << 62;
+ else
+ bit = 0x0000004000000000LL; //(sll)1 << 38;
+
+ while (bit > num) bit >>= 2;
+
+ while (bit)
+ {
+ if (num >= result + bit)
+ {
+ num -= result + bit;
+ result = (result >> 1) + bit;
+ }
+ else
+ {
+ result = (result >> 1);
+ }
+ bit >>= 2;
+ }
+
+ // Then process it again to get the lowest 8 bits.
+ // 尾数定点数大于等于1
+ if (num > CONST_P9999)
+ {
+ // The remainder 'num' is too large to be shifted left
+ // by 32, so we have to add 0.5 to result manually and
+ // adjust 'num' accordingly. a是原本的数
+ // why 0.5 is enough? because: result^2 <= a < (result+1)^2
+ // num = a - (result + 0.5)^2
+ // = num + result^2 - (result + 0.5)^2
+ // = num - result - 0.25
+ num -= result;
+ // 注意，严格来说，这里要左移16位才缩小到结果，但是为了后续小数不用再放大，这里提前做了放大
+ num = (num << 32) - CONST_1_4;
+ result = (result << 32) + CONST_1_2;
+ }
+ else
+ {
+ num <<= 32;
+ result <<= 32;
+ }
+
+ bit = 0x0000000040000000LL; //(sll)1 << 30;
+
+ while (bit)
+ {
+ if (num >= result + bit)
+ {
+ num -= result + bit;
+ result = (result >> 1) + bit;
+ }
+ else
+ {
+ result = (result >> 1);
+ }
+ bit >>= 2;
+ }
+
+ return result;
+} \ No newline at end of file
diff --git a/Client/ThirdParty/math-sll/math-sll.h b/Client/ThirdParty/math-sll/math-sll.h
new file mode 100644
index 0000000..b10dbe7
--- /dev/null
+++ b/Client/ThirdParty/math-sll/math-sll.h
@@ -0,0 +1,770 @@
+#if !defined(MATH_SLL_H)
+# define MATH_SLL_H
+
+/*
+ * Revision v1.24
+ *
+ * A fixed point (32.32 bit) math library.
+ *
+ * Description
+ *
+ * Floating point packs the most accuracy in the available bits, but it
+ * often provides more accuracy than is required. It is time consuming to
+ * carry the extra precision around, particularly on platforms that don't
+ * have a dedicated floating point processor.
+ *
+ * This library is a compromise. All math is done using the 64 bit signed
+ * "long long" format (sll), and is not intended to be portable, just as
+ * simple and as fast as possible.
+ *
+ * As some processors lack division instructions but have multiplication
+ * instructions, multiplication is favored over division. This can be a
+ * penalty when used on a processor with a division instruction, so it is
+ * recommended to modify the division functions and macros in that case.
+ *
+ * On procesors without multiplication instructions, other algorithms, for
+ * example CORDIC, are probably faster.
+ *
+ * Since "long long" is a elementary type, it can be passed around without
+ * resorting to the use of pointers. Since the format used is fixed point,
+ * there is never a need to do time consuming checks and adjustments to
+ * maintain normalized numbers, as is the case in floating point.
+ *
+ * Simply put, this library is limited to handling numbers with a whole
+ * part of up to 2^31 - 1 = 2.147483647e9 in magnitude, and fractional
+ * parts down to 2^-32 = 2.3283064365e-10 in magnitude. This yields a
+ * decent range and accuracy for many applications.
+ *
+ * IMPORTANT
+ *
+ * No checking for arguments out of range (error).
+ * No checking for divide by zero (error).
+ * No checking for overflow (error).
+ * No checking for underflow (warning).
+ * Chops, doesn't round.
+ *
+ * Functions
+ *
+ * sll dbl2sll(double d) double to sll
+ * double sll2dbl(sll s) sll to double
+ *
+ * sll int2sll(int i) integer to sll
+ * int sll2int(sll s) sll to integer
+ *
+ * sll sllint(sll s) Set fractional-part to 0
+ * sll sllfrac(sll s) Set integer-part to 0
+ *
+ * sll slladd(sll x, sll y) x + y
+ * sll sllneg(sll x) -x
+ * sll sllsub(sll x, sll y) x - y
+ *
+ * sll sllmul(sll x, sll y) x * y
+ * sll sllmul2(sll x) x * 2
+ * sll sllmul2n(sll x, int n) x * 2^n, 0 <= n <= 31
+ * sll sllmul4(sll x) x * 4
+ *
+ * sll slldiv(sll x, sll y) x / y
+ * sll slldiv2(sll x) x / 2
+ * sll slldiv2n(sll x, int n) x / 2^n, 0 <= n <= 31
+ * sll slldiv4(sll x) x / 4
+ *
+ * sll sllcos(sll x) cos x
+ * sll sllsin(sll x) sin x
+ * sll slltan(sll x) tan x
+ *
+ * sll sllsec(sll x) sec x = 1 / cos x
+ * sll sllcsc(sll x) csc x = 1 / sin x
+ * sll sllcot(sll x) cot x = 1 / tan x = cos x / sin x
+ *
+ * sll sllacos(sll x) acos x
+ * sll sllasin(sll x) asin x
+ * sll sllatan(sll x) atan x
+ *
+ * sll sllcosh(sll x) cosh x
+ * sll sllsinh(sll x) sinh x
+ * sll slltanh(sll x) tanh x
+ *
+ * sll sllsech(sll x) sech x
+ * sll sllcsch(sll x) cosh x
+ * sll sllcoth(sll x) coth x
+ *
+ * sll sllexp(sll x) e^x
+ * sll slllog(sll x) ln x
+ *
+ * sll sllinv(sll v) 1 / x
+ * sll sllpow(sll x, sll y) x^y
+ * sll sllsqrt(sll x) x^(1 / 2)
+ *
+ * sll sllfloor(sll x) floor x
+ * sll sllceil(sll x) ceiling x
+ *
+ * Macros
+ *
+ * Use of the following macros is optional, but may be beneficial with
+ * some compilers. Using the non-macro versions is strongly recommended.
+ *
+ * WARNING: macros do not type-check their arguments!
+ *
+ * _int2sll(X) See function int2sll()
+ * _sll2int(X) See function sll2int()
+ *
+ * _sllint(X) See function sllint()
+ * _sllfrac(X) See function sllfrac()
+ *
+ * _slladd(X,Y) See function slladd()
+ * _sllneg(X) See function sllneg()
+ * _sllsub(X,Y) See function sllsub()
+ *
+ * _sllmul2(X) See function sllmul2()
+ * _sllmul2n(X) See function sllmul2n()
+ * _sllmul4(X) See function sllmul4()
+ *
+ * _slldiv(X,Y) See function slldiv()
+ * _slldiv2(X,Y) See function slldiv2()
+ * _slldiv2n(X,Y) See function slldiv2n()
+ * _slldiv4(X,Y) See function slldiv4()
+ *
+ * Credits
+ *
+ * Maintained, conceived, written, and fiddled with by:
+ *
+ * Andrew E. Mileski <andrewm@isoar.ca>
+ *
+ * Other source code contributors:
+ *
+ * Kevin Rockel
+ * Kevin Michael Woley
+ * Mark Anthony Lisher
+ * Nicolas Pitre
+ * Anonymous
+ *
+ * License
+ *
+ * Licensed under the terms of the MIT license:
+ *
+ * Copyright (c) 2000,2002,2006,2012,2016 Andrew E. Mileski <andrewm@isoar.ca>
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to
+ * deal in the Software without restriction, including without limitation the
+ * rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
+ * sell copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The copyright notice, and this permission notice shall be included in all
+ * copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+ * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+ * DEALINGS IN THE SOFTWARE.
+ */
+
+/* DEC SA-110 "StrongARM" (armv4l) architecture has a big-endian double */
+#if defined(__arm__)
+# if (!defined(__BYTE_ORDER__) || (__BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__))
+# define BROKEN_IEEE754_DOUBLE
+# endif
+#endif
+
+#ifdef _MSC_VER
+# define __inline__ __inline
+# ifndef _MSC_STDINT_H_
+typedef signed __int64 int64_t;
+typedef unsigned __int64 uint64_t;
+#endif
+#define __extension__
+#else
+#include <stdint.h>
+#if defined(__GNUC__)
+
+#else
+ #define __inline__ inline
+#endif
+#endif
+/*
+ * Data types
+ */
+
+__extension__ typedef int64_t sll;
+__extension__ typedef uint64_t ull;
+
+/*
+ * Function prototypes
+ */
+
+sll dbl2sll(double d);
+double sll2dbl(sll s);
+
+static __inline__ sll int2sll(int i);
+static __inline__ int sll2int(sll s);
+
+static __inline__ sll sllint(sll s);
+static __inline__ sll sllfrac(sll s);
+
+static __inline__ sll slladd(sll x, sll y);
+static __inline__ sll sllneg(sll s);
+static __inline__ sll sllsub(sll x, sll y);
+sll sllmul(sll x, sll y);
+static __inline__ sll sllmul2(sll x);
+static __inline__ sll sllmul4(sll x);
+static __inline__ sll sllmul2n(sll x, int n);
+
+static __inline__ sll slldiv(sll x, sll y);
+static __inline__ sll slldiv2(sll x);
+static __inline__ sll slldiv4(sll x);
+static __inline__ sll slldiv2n(sll x, int n);
+
+sll sllcos(sll x);
+sll sllsin(sll x);
+sll slltan(sll x);
+
+static __inline__ sll sllacos(sll x);
+sll sllasin(sll x);
+sll sllatan(sll x);
+
+static __inline__ sll sllsec(sll x);
+static __inline__ sll sllcsc(sll x);
+static __inline__ sll sllcot(sll x);
+
+static __inline__ sll sllcosh(sll x);
+static __inline__ sll sllsinh(sll x);
+static __inline__ sll slltanh(sll x);
+
+static __inline__ sll sllsech(sll x);
+static __inline__ sll sllcsch(sll x);
+static __inline__ sll sllcoth(sll x);
+
+sll sllexp(sll x);
+sll slllog(sll x);
+
+sll sllpow(sll x, sll y);
+sll sllinv(sll v);
+sll sllsqrt(sll x);
+sll slld2dsqrt(sll x);
+
+static __inline__ sll sllfloor(sll x);
+static __inline__ sll sllceil(sll x);
+
+/*
+ * Macros
+ *
+ * WARNING - Macros don't type-check!
+ */
+
+#define _int2sll(X) (((sll) (X)) << 32)
+#define _sll2int(X) ((int) ((X) >> 32))
+
+#define _sllint(X) ((X) & 0xffffffff00000000LL)
+#define _sllfrac(X) ((X) & 0x00000000ffffffffLL)
+
+#define _slladd(X,Y) ((X) + (Y))
+#define _sllneg(X) (-(X))
+#define _sllsub(X,Y) ((X) - (Y))
+
+#define _sllmul2(X) ((X) << 1)
+#define _sllmul4(X) ((X) << 2)
+#define _sllmul2n(X,N) ((X) << (N))
+
+#define _slldiv(X,Y) sllmul((X), sllinv(Y))
+#define _slldiv2(X) ((X) >> 1)
+#define _slldiv4(X) ((X) >> 2)
+#define _slldiv2n(X,N) ((X) >> (N))
+
+/*
+ * Constants (converted from double)
+*/
+
+#define CONST_0 0x0000000000000000LL // 0.0
+#define CONST_1 0x0000000100000000LL // 1.0
+#define CONST_neg1 0xffffffff00000000LL // -1.0
+#define CONST_2 0x0000000200000000LL // 2.0
+#define CONST_3 0x0000000300000000LL // 3.0
+#define CONST_4 0x0000000400000000LL // 4.0
+#define CONST_10 0x0000000a00000000LL // 10.0
+#define CONST_10 0x0000000a00000000LL // 10.0
+#define CONST_1_2 0x0000000080000000LL // 1.0 / 2.0
+#define CONST_1_3 0x0000000055555555LL // 1.0 / 3.0
+#define CONST_1_4 0x0000000040000000LL // 1.0 / 4.0
+#define CONST_1_5 0x0000000033333333LL // 1.0 / 5.0
+#define CONST_1_6 0x000000002aaaaaaaLL // 1.0 / 6.0
+#define CONST_1_7 0x0000000024924924LL // 1.0 / 7.0
+#define CONST_1_8 0x0000000020000000LL // 1.0 / 8.0
+#define CONST_1_9 0x000000001c71c71cLL // 1.0 / 9.0
+#define CONST_1_10 0x0000000019999999LL // 1.0 / 10.0
+#define CONST_1_11 0x000000001745d174LL // 1.0 / 11.0
+#define CONST_1_12 0x0000000015555555LL // 1.0 / 12.0
+#define CONST_1_20 0x000000000cccccccLL // 1.0 / 20.0
+#define CONST_1_30 0x0000000008888888LL // 1.0 / 30.0
+#define CONST_1_42 0x0000000006186186LL // 1.0 / 42.0
+#define CONST_1_56 0x0000000004924924LL // 1.0 / 56.0
+#define CONST_1_72 0x00000000038e38e3LL // 1.0 / 72.0
+#define CONST_1_90 0x0000000002d82d82LL // 1.0 / 90.0
+#define CONST_1_110 0x000000000253c825LL // 1.0 / 110.0
+#define CONST_1_132 0x0000000001f07c1fLL // 1.0 / 132.0
+#define CONST_1_156 0x0000000001a41a41LL // 1.0 / 156.0
+#define CONST_P9999 0x00000000FFFFFFFFLL // 0.999999999999
+
+#define CONST_E 0x00000002b7e15162LL // E
+#define CONST_1_E 0x000000005e2d58d8LL // 1 / E
+#define CONST_SQRTE 0x00000001a61298e1LL // sqrt(E)
+#define CONST_1_SQRTE 0x000000009b4597e3LL // 1 / sqrt(E)
+#define CONST_LOG2_E 0x0000000171547652LL // ln(E)
+#define CONST_LOG10_E 0x000000006f2dec54LL // log(E)
+#define CONST_LN2 0x00000000b17217f7LL // ln(2)
+#define CONST_LN10 0x000000024d763776LL // ln(10)
+
+#define CONST_PI 0x00000003243f6a88LL // PI
+#define CONST_2PI 0x00000006487ED510LL // PI
+#define CONST_PI_2 0x00000001921fb544LL // PI / 2
+#define CONST_PI_4 0x00000000c90fdaa2LL // PI / 4
+#define CONST_1_PI 0x00000000517cc1b7LL // 1 / PI
+#define CONST_2_PI 0x00000000a2f9836eLL // 2 / PI
+#define CONST_180_PI 0x000000394BB834C7LL // 180 / PI 246083499207.51537232162612011973
+#define CONST_PI_180 0x000000000477d1a8LL // PI / 180 74961320.580677883103327681382757
+#define CONST_2_SQRTPI 0x0000000120dd7504LL // 2 / sqrt(PI)
+#define CONST_SQRT2 0x000000016a09e667LL // sqrt(2)
+#define CONST_1_SQRT2 0x00000000b504f333LL // 1 / sqrt(2)
+
+#define CONST_FACT_0 0x0000000100000000LL // 0!
+#define CONST_FACT_1 0x0000000100000000LL // 1!
+#define CONST_FACT_2 0x0000000200000000LL // 2!
+#define CONST_FACT_3 0x0000000600000000LL // 3!
+#define CONST_FACT_4 0x0000001800000000LL // 4!
+#define CONST_FACT_5 0x0000007800000000LL // 5!
+#define CONST_FACT_6 0x000002d000000000LL // 6!
+#define CONST_FACT_7 0x000013b000000000LL // 7!
+#define CONST_FACT_8 0x00009d8000000000LL // 8!
+#define CONST_FACT_9 0x0005898000000000LL // 9!
+#define CONST_FACT_10 0x00375f0000000000LL // 10!
+#define CONST_FACT_11 0x0261150000000000LL // 11!
+#define CONST_FACT_12 0x1c8cfc0000000000LL // 12!
+
+#define CONST_MAX 0x7FFFFFFFFFFFFFFFLL // 最大
+#define CONST_MIN 0x8000000000000000LL // 最小
+
+/*
+ * Convert integer to sll
+ */
+
+static __inline__ sll int2sll(int i)
+{
+ return _int2sll(i);
+}
+
+/*
+ * Convert sll to integer (truncates)
+ */
+
+static __inline__ int sll2int(sll s)
+{
+ return _sll2int(s);
+}
+
+/*
+ * Integer-part of sll (fractional-part set to 0)
+ */
+
+static __inline__ sll sllint(sll s)
+{
+ return _sllint(s);
+}
+
+/*
+ * Fractional-part of sll (integer-part set to 0)
+ */
+
+static __inline__ sll sllfrac(sll s)
+{
+ return _sllfrac(s);
+}
+
+/*
+ * Addition
+ */
+
+static __inline__ sll slladd(sll x, sll y)
+{
+ return _slladd(x, y);
+}
+
+/*
+ * Negate
+ */
+
+static __inline__ sll sllneg(sll s)
+{
+ return _sllneg(s);
+}
+
+/*
+ * Subtraction
+ */
+
+static __inline__ sll sllsub(sll x, sll y)
+{
+ return _sllsub(x, y);
+}
+
+/*
+ * Multiply two sll values
+ *
+ * Description
+ *
+ * Let a = A * 2^32 + a_h * 2^0 + a_l * 2^(-32)
+ * Let b = B * 2^32 + b_h * 2^0 + b_l * 2^(-32)
+ *
+ * Where:
+ *
+ * *_h is the integer part
+ * *_l the fractional part
+ * A and B are the sign (0 for positive, -1 for negative).
+ *
+ * a * b = (A * 2^32 + a_h * 2^0 + a_l * 2^-32)
+ * * (B * 2^32 + b_h * 2^0 + b_l * 2^-32)
+ *
+ * Expanding the terms, we get:
+ *
+ * = A * B * 2^64 + A * b_h * 2^32 + A * b_l * 2^0
+ * + a_h * B * 2^32 + a_h * b_h * 2^0 + a_h * b_l * 2^-32
+ * + a_l * B * 2^0 + a_l * b_h * 2^-32 + a_l * b_l * 2^-64
+ *
+ * Grouping by powers of 2, we get:
+ *
+ * = A * B * 2^64
+ * Meaningless overflow from sign extension - ignore
+ *
+ * + (A * b_h + a_h * B) * 2^32
+ * Overflow which we can't handle - ignore
+ *
+ * + (A * b_l + a_h * b_h + a_l * B) * 2^0
+ * We only need the low 32 bits of this term, as the rest is overflow
+ *
+ * + (a_h * b_l + a_l * b_h) * 2^-32
+ * We need all 64 bits of this term
+ *
+ * + a_l * b_l * 2^-64
+ * We only need the high 32 bits of this term, as the rest is underflow
+ *
+ * Note that:
+ * a > 0 && b > 0: A = 0, B = 0 and the third term is a_h * b_h
+ * a < 0 && b > 0: A = -1, B = 0 and the third term is a_h * b_h - b_l
+ * a > 0 && b < 0: A = 0, B = -1 and the third term is a_h * b_h - a_l
+ * a < 0 && b < 0: A = -1, B = -1 and the third term is a_h * b_h - a_l - b_l
+ */
+
+/*
+ * Multiplication by 2
+ */
+
+static __inline__ sll sllmul2(sll x)
+{
+ return _sllmul2(x);
+}
+
+/*
+ * Multiplication by 4
+ */
+
+static __inline__ sll sllmul4(sll x)
+{
+ return _sllmul4(x);
+}
+
+/*
+ * Multiplication by power of 2
+ */
+
+static __inline__ sll sllmul2n(sll x, int n)
+{
+ return _sllmul2n(x, n);
+}
+
+/*
+ * Division
+ */
+
+static __inline__ sll slldiv(sll x, sll y)
+{
+ return _slldiv(x, y);
+}
+
+/*
+ * Division by 2
+ */
+
+static __inline__ sll slldiv2(sll x)
+{
+ return _slldiv2(x);
+}
+
+/*
+ * Division by 4
+ */
+
+static __inline__ sll slldiv4(sll x)
+{
+ return _slldiv4(x);
+}
+
+/*
+ * Division by power of 2
+ */
+
+static __inline__ sll slldiv2n(sll x, int n)
+{
+ return _slldiv2n(x, n);
+}
+
+/*
+ *
+ * Calculate acos x, where |x| <= 1
+ *
+ * Description
+ *
+ * acos x = pi / 2 - asin x
+ * acos x = pi / 2 - SUM[n=0,) C(2 * n, n) * x^(2 * n + 1) / (4^n * (2 * n + 1)), |x| <= 1
+ *
+ * where C(n, r) = nCr = n! / (r! * (n - r)!)
+ */
+
+static __inline__ sll sllacos(sll x)
+{
+ return _sllsub(CONST_PI_2, sllasin(x));
+}
+
+/*
+ * Trigonometric secant
+ *
+ * Description
+ *
+ * sec x = 1 / cos x
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllsec(sll x)
+{
+ return sllinv(sllcos(x));
+}
+
+/*
+ * Trigonometric cosecant
+ *
+ * Description
+ *
+ * csc x = 1 / sin x
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllcsc(sll x)
+{
+ return sllinv(sllsin(x));
+}
+
+/*
+ * Trigonometric cotangent
+ *
+ * Description
+ *
+ * cot x = 1 / tan x
+ *
+ * cot x = cos x / sin x
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllcot(sll x)
+{
+ return _slldiv(sllcos(x), sllsin(x));
+}
+
+/*
+ * Hyperbolic cosine
+ *
+ * Description
+ *
+ * cosh x = (e^x + e^(-x)) / 2
+ *
+ * cosh x = 1 + x^2 / 2! + ... + x^(2 * N) / (2 * N)!
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllcosh(sll x)
+{
+ return _slldiv2(_slladd(sllexp(x), sllexp(_sllneg(x))));
+}
+
+/*
+ * Hyperbolic sine
+ *
+ * Description
+ *
+ * sinh x = (e^x - e^(-x)) / 2
+ *
+ * sinh x = 1 + x^3 / 3! + ... + x^(2 * N + 1) / (2 * N + 1)!
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllsinh(sll x)
+{
+ return _slldiv2(_sllsub(sllexp(x), sllexp(_sllneg(x))));
+}
+
+/*
+ * Hyperbolic tangent
+ *
+ * Description
+ *
+ * tanh x = sinh x / cosh x
+ *
+ * tanh x = (e^x - e^(-x)) / (e^x + e^(-x))
+ *
+ * tanh x = (e^(2 * x) - 1) / (e^(2 * x) + 1)
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll slltanh(sll x)
+{
+ register sll e2x;
+
+ e2x = sllexp(_sllmul2(x));
+
+ return _slldiv(_sllsub(e2x, CONST_1), _slladd(e2x, CONST_1));
+}
+
+/*
+ * Hyperbolic secant
+ *
+ * Description
+ *
+ * sech x = 1 / cosh x
+ *
+ * sech x = 2 / (e^x + e^(-x))
+ *
+ * sech x = 2 * e^x / (e^(2 * x) + 1)
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllsech(sll x)
+{
+ return _slldiv(_sllmul2(sllexp(x)), _slladd(sllexp(_sllmul2(x)), CONST_1));
+}
+
+/*
+ * Hyperbolic cosecant
+ *
+ * Description
+ *
+ * csch x = = 1 / sinh x
+ *
+ * csch x = 2 / (e^x - e^(-x))
+ *
+ * csch x = 2 * e^x / (e^(2 * x) - 1)
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllcsch(sll x)
+{
+ return _slldiv(_sllmul2(sllexp(x)), _sllsub(sllexp(_sllmul2(x)), CONST_1));
+}
+
+/*
+ * Hyperbolic cotangent
+ *
+ * Description
+ *
+ * coth x = 1 / tanh x
+ *
+ * coth x = cosh x / sinh x
+ *
+ * coth x = (e^x + e^(-x)) / (e^x - e^(-x))
+ *
+ * coth x = (e^(2 * x) + 1) / (e^(2 * x) - 1)
+ *
+ * An alternate algorithm, like a power series, would be more accurate.
+ */
+
+static __inline__ sll sllcoth(sll x)
+{
+ register sll e2x;
+
+ e2x = sllexp(sllmul2(x));
+
+ return _slldiv(_slladd(e2x, CONST_1), _sllsub(e2x, CONST_1));
+}
+
+/*
+ * Floor
+ *
+ * Description
+ *
+ * floor x = largest integer not larger than x
+ */
+
+static __inline__ sll sllfloor(sll x)
+{
+ register sll retval;
+
+ retval = _sllint(x);
+
+ return ((retval > x) ? _sllsub(retval, CONST_1): retval);
+}
+
+/*
+ * Ceiling
+ *
+ * Description
+ *
+ * ceil x = smallest integer not smaller than x
+ */
+
+static __inline__ sll sllceil(sll x)
+{
+ register sll retval;
+
+ retval = _sllint(x);
+
+ return ((retval < x) ? _slladd(retval, CONST_1): retval);
+}
+
+#define sllabs(x) ((x) < 0 ? -(x) : (x))
+#define __VECTOR2_META__ "__VECTOR2_META__"
+#define __VECTOR3_META__ "__VECTOR3_META__"
+#define __ROT2_META__ "__ROT2_META__"
+#define __ROT4_META__ "__ROT4_META__"
+#define __METATABLE_NAME "__FIX_METATABLE__"
+
+static int _mul[] = {1, 10, 100, 1000, 10000, 100000, 1000000};
+
+typedef struct Vector2
+{
+ sll x;
+ sll y;
+}Vector2;
+
+typedef struct Vector3
+{
+ sll x;
+ sll y;
+ sll z;
+}Vector3;
+
+typedef struct Vector4
+{
+ sll x;
+ sll y;
+ sll z;
+ sll w;
+
+}Vector4;
+#endif
generated by cgit v1.1-26-g67d0 at 2025-05-29 09:00:31 +0000