diff --git a/Tests/Cosmos.Compiler.Tests.Bcl/System/DoubleTest.cs b/Tests/Cosmos.Compiler.Tests.Bcl/System/DoubleTest.cs index ed8dd6e50..5188452b1 100644 --- a/Tests/Cosmos.Compiler.Tests.Bcl/System/DoubleTest.cs +++ b/Tests/Cosmos.Compiler.Tests.Bcl/System/DoubleTest.cs @@ -187,6 +187,31 @@ namespace Cosmos.Compiler.Tests.Bcl.System value = 42.0; valueNegated = -value; Assert.IsTrue((EqualityHelper.DoublesAreEqual(valueNegated, -42.0f)), "(double) negation of positive double doesn't work"); + + #region Parsing + value = double.Parse("0.4"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.4), "simple parsing of double works"); + + value = double.Parse("+0.3"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.3), "parsing of double with positive sign works!"); + + value = double.Parse("-0.4"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, -0.4), "parsing of negative double works!"); + + value = double.Parse(" 0.7 "); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.7), "double parsing ignores leading and trailing whitespaces"); + + value = double.Parse("0.4E1"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 4), "double parsing takes in account E"); + + value = double.Parse("0.4E-1"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.04), "double parsing works with negative E"); + + Assert.IsFalse(double.TryParse("asd4", out value), "double TryParse returns false when it fails"); + + Assert.IsTrue(double.TryParse("2.3", out value), "double TryParse returns true when it works"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 2.3), "double TryParse returns correct result when it works"); + #endregion } } } diff --git a/Tests/Cosmos.Compiler.Tests.Bcl/System/MathTest.cs b/Tests/Cosmos.Compiler.Tests.Bcl/System/MathTest.cs index eb02c05d5..ddebd21c3 100644 --- a/Tests/Cosmos.Compiler.Tests.Bcl/System/MathTest.cs +++ b/Tests/Cosmos.Compiler.Tests.Bcl/System/MathTest.cs @@ -6,7 +6,7 @@ using Cosmos.Compiler.Tests.Bcl.Helper; namespace Cosmos.Compiler.Tests.Bcl.System { - class MathTest + internal class MathTest { public static void Execute() { @@ -35,6 +35,126 @@ namespace Cosmos.Compiler.Tests.Bcl.System // Test with positive infinity result = Math.Sqrt(double.PositiveInfinity); Assert.IsTrue(double.IsPositiveInfinity(result), "Sqrt of PositiveInfinity must return PositiveInfinity"); + + #region Math.Exp + + //Test with integer + result = Math.Exp(2); + Assert.IsTrue((result == 7.38905609893065), "e^2 is equal to 7.38905609893065"); + + //Test with double exponent + result = Math.Exp(1.5); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 4.48168907033806), "e^1.5 returns correct result"); + + result = Math.Exp(0); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 1), "e^0 gives correct result"); + + result = Math.Exp(1); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, Math.E), "e^1 gives correct result"); + + result = Math.Exp(double.PositiveInfinity); + Assert.IsTrue(result == double.PositiveInfinity, "e^Infinity gives correct result"); + + result = Math.Exp(double.NegativeInfinity); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "e^-Infinity gives correct result"); + + result = Math.Exp(double.NaN); + Assert.IsTrue(double.IsNaN(result), "e^NaN gives correct result"); + + result = Math.Exp(double.MaxValue); + Assert.IsTrue(double.IsPositiveInfinity(result), "e^0 gives correct result"); + + result = Math.Exp(double.MinValue); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "e^0 gives correct result"); + + #endregion Math.Exp + + #region Math.Pow + + //Test with integer power + result = Math.Pow(2, 2); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 4), "2^2 gives accurate result"); + + //Test with decimal power + result = Math.Pow(9, 0.5); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, Math.Sqrt(9)), "9^0.5 gives same answer as sqrt(9)"); + + //Test with negative base + result = Math.Pow(-2, 2); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 4), "Math.Pow gives correct result when raising negative number to even power"); + + result = Math.Pow(-2, 3); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, -8), "Math.Pow gives correct result when raising negative number to odd power"); + + //Test with negative power + result = Math.Pow(2, -1); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0.5), "Pow gives correct results when handling negative powers"); + + //Have double as base + result = Math.Pow(0.5, 2); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0.25), "Pow gives correct solution with double base"); + + //x = Nan + result = Math.Pow(double.NaN, 2); + Assert.IsTrue(double.IsNaN(result), "Pow gives correct solution when x is NaN"); + //Y = Nan + result = Math.Pow(10, double.NaN); + Assert.IsTrue(double.IsNaN(result), "Pow gives correct solution when y is NaN"); + //y = 0 + result = Math.Pow(10, 0); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 1), "Pow gives correct solution when y is 0"); + //x = -Inf y < 0 == 0 + result = Math.Pow(double.NegativeInfinity, -2); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "Pow gives correct solution when X is -INF and y is negative"); + //x = -Inf y > 0 && y is even == Inf + result = Math.Pow(double.NegativeInfinity, 2); + Assert.IsTrue(double.IsPositiveInfinity(result), "Pow gives correct solution when x is -INF and y is even"); + //x is -INF and y is positive odd == -INF + result = Math.Pow(double.NegativeInfinity, 3); + Assert.IsTrue(double.IsNegativeInfinity(result), "Pow gives correct solution when x is -INF and y is odd"); + //x < 0 && y is not integer or special case + result = Math.Pow(-3, 0.25); + Assert.IsTrue(double.IsNaN(result), "Pow gives correct solution when x is negative and y is non integer"); + //x = -1 && y is Inf == Nan + result = Math.Pow(-1, double.PositiveInfinity); + Assert.IsTrue(double.IsNaN(result), "Pow gives correct solution when x is -1 and y is INF"); + //x = -1 && y is -Inf == Nan + result = Math.Pow(-1, double.NegativeInfinity); + Assert.IsTrue(double.IsNaN(result), "Pow gives correct solution when x is -1 and y is -INF"); + //-1 < x < 1 + y = -Inf + result = Math.Pow(-0.25, double.NegativeInfinity); + Assert.IsTrue(double.IsPositiveInfinity(result), "Pow gives correct solution when -1 < x < 0 and y = -INF"); + result = Math.Pow(0.25, double.NegativeInfinity); + Assert.IsTrue(double.IsPositiveInfinity(result), "Pow gives correct solution when 0 < x < 1 and y = -INF"); + //-1 < x < 1 + y = Inf + result = Math.Pow(-0.25, double.PositiveInfinity); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "Pow gives correct solution when -1 < x < 0 and y is INF"); + result = Math.Pow(0.25, double.PositiveInfinity); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "Pow gives correct solution when 0 < x < 1 and y is INF"); + //-1 > x || x > 1 + y = -Inf + result = Math.Pow(-1.5, double.NegativeInfinity); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "Pow gives correct solution when x < -1 and y is -INF"); + result = Math.Pow(1.5, double.NegativeInfinity); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "Pow gives correct solution when x > 1 and y is -INF"); + //-1 > x || x > 1 + y = Inf + result = Math.Pow(-1.25, double.PositiveInfinity); + Assert.IsTrue(double.IsPositiveInfinity(result), "Pow gives correct solution when -1 > x and y = INF"); + result = Math.Pow(1.25, double.PositiveInfinity); + Assert.IsTrue(double.IsPositiveInfinity(result), "Pow gives correct solution when x > 1 and y = INF"); + //x = 0 y > 0 + result = Math.Pow(0, 2); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "Pow gives correct solution when x = 0 any y > 0 "); + //x = 0 y < 0 + result = Math.Pow(0, -3); + Assert.IsTrue(double.IsPositiveInfinity(result), "Pow gives correct solution when x is 0 and y < 0 "); + //x = inf y < 0 + result = Math.Pow(double.PositiveInfinity, -5); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(result, 0), "Pow gives correct solution when x is INF and y < 0 "); + //x = inf y > 0 + result = Math.Pow(double.PositiveInfinity, 5); + Assert.IsTrue(double.IsPositiveInfinity(result), "Pow gives correct solution when x is INF and y > 0 "); + + #endregion Math.Pow } } } diff --git a/Tests/Cosmos.Compiler.Tests.Bcl/System/SingleTest.cs b/Tests/Cosmos.Compiler.Tests.Bcl/System/SingleTest.cs index 9f86a9a44..3efdb1e27 100644 --- a/Tests/Cosmos.Compiler.Tests.Bcl/System/SingleTest.cs +++ b/Tests/Cosmos.Compiler.Tests.Bcl/System/SingleTest.cs @@ -1,17 +1,19 @@ using System; using System.Linq; using System.Threading.Tasks; +using Cosmos.Compiler.Tests.Bcl.Helper; using Cosmos.TestRunner; namespace Cosmos.Compiler.Tests.Bcl.System { - class SingleTest + internal class SingleTest { /* The single== equality operator is so imprecise to not be really ever useful we should be happy if the two values are "similar" */ + private static bool SinglesAreEqual(Single left, Single right) { // Define the tolerance for variation in their values - Single difference = (Single) Math.Abs(left * .00001); + Single difference = (Single)Math.Abs(left * .00001); if (Math.Abs(left - right) <= difference) return true; @@ -62,7 +64,7 @@ namespace Cosmos.Compiler.Tests.Bcl.System expectedResult = "-42.42"; Assert.IsTrue((result == expectedResult), "Single.ToString of negative number doesn't work"); - + /* A big value (to be correct toString should convert it in scientific notation) */ value = 9223372036854775808f; @@ -73,7 +75,7 @@ namespace Cosmos.Compiler.Tests.Bcl.System /* OK now a normal value */ value = 42.42F; // It exists Single.MaxValue but it is a too big value an can be represented only on Scientific notation but then how to confront with a String? - + result = value.ToString(); expectedResult = "42.42"; @@ -206,6 +208,33 @@ namespace Cosmos.Compiler.Tests.Bcl.System value = 42.0f; valueNegated = -value; Assert.IsTrue((SinglesAreEqual(valueNegated, -42.0f)), "(float) negation of positive float doesn't work"); + + #region Parsing + + value = float.Parse("0.4"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.4), "simple parsing of float works"); + + value = float.Parse("+0.3"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.3), "parsing of float with positive sign works!"); + + value = float.Parse("-0.4"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, -0.4), "parsing of negative float works!"); + + value = float.Parse(" 0.7 "); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.7), "float parsing ignores leading and trailing whitespaces"); + + value = float.Parse("0.4E1"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 4), "float parsing takes in account E"); + + value = float.Parse("0.4E-1"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 0.04), "float parsing works with negative E"); + + Assert.IsFalse(float.TryParse("asd4", out value), "float TryParse returns false when it fails"); + + Assert.IsTrue(float.TryParse("2.3", out value), " float TryParse returns true when it works"); + Assert.IsTrue(EqualityHelper.DoublesAreEqual(value, 2.3), "float TryParse returns correct result when it works"); + + #endregion Parsing } } } diff --git a/source/Cosmos.System2_Plugs/System/DoubleImpl.cs b/source/Cosmos.System2_Plugs/System/DoubleImpl.cs index 881cc2c6c..c9457b216 100644 --- a/source/Cosmos.System2_Plugs/System/DoubleImpl.cs +++ b/source/Cosmos.System2_Plugs/System/DoubleImpl.cs @@ -1,16 +1,118 @@ using System; using Cosmos.Common; -using IL2CPU.API; +using System.Collections.Generic; using IL2CPU.API.Attribs; namespace Cosmos.System_Plugs.System { - [Plug(Target = typeof (double))] + [Plug(Target = typeof(double))] public static class DoubleImpl { public static string ToString(ref double aThis) { return StringHelper.GetNumberString(aThis); } + + public static double Parse(string s) + { + //Format of Double string: [whitespace][sign][integral-digits[,]]integral-digits[.[fractional-digits]][E[sign]exponential-digits][whitespace] + + //Validate input + if (s is null) throw new ArgumentNullException("s can not be null"); + + //Remove leading whitespaces + while (s.Length != 0 && (s[0] == ' ' || s[0] == '\n' || s[0] == '\t')) + { + s = s.Substring(1); + } + + //Check that string is not finished too early + if (s.Length == 0) throw new FormatException(); + + //Check for sign + short sign = 1; + if (s[0] == '-') + { + s = s.Substring(1); + sign = -1; + } + else if (s[0] == '+') s = s.Substring(1); + + //Check that string is not finished too early + if (s.Length == 0) throw new FormatException(); + + //Read in number + + List internalDigits = new List(); + List fractionalDigits = new List(); + + bool foundDecimal = false; + + //Iterate until fully parsed or an E/Whitespace is found + //Assume english standard, so . == decimal seperator and , == thousands seperator + while (s.Length != 0) + { + char active = s[0]; + if (active == 'E' || active == 'e' || active == ' ' || active == '\n' || active == '\t') break; + + s = s.Substring(1); + if (active == '.') foundDecimal = true; + else if (active == ',') continue; + else if (active >= '0' && active <= '9') + { + if (foundDecimal) fractionalDigits.Add(int.Parse(active.ToString())); + else internalDigits.Add(int.Parse(active.ToString())); + } + else + { + throw new FormatException(); + } + } + + //Iterate through rest of string + double multiplier = 0; + while (s.Length != 0) + { + //Check for exponential notation i.e. 8.1E10 = 8.1 * 10^10 + Whitespaces + //E can only be followed by integers + if (s[0] == 'E' || s[0] == 'e') + { + multiplier = double.Parse(s.Substring(1)); + break; + } + else if (s[0] == ' ' || s[0] == '\n' || s[0] == '\t') s = s.Substring(1); + else throw new FormatException(); + } + + //Create double + double parsed = 0; + for (int i = 0; i < internalDigits.Count; i++) + { + parsed += internalDigits[i] * Math.Pow(10, (internalDigits.Count - (i + 1))); + } + for (int i = 0; i < fractionalDigits.Count; i++) + { + parsed += fractionalDigits[i] * Math.Pow(10, -1 * (i + 1)); + } + + parsed *= Math.Pow(10, multiplier); + parsed *= sign; + + return parsed; + } + + public static bool TryParse(string s, out double result) + { + try + { + result = Parse(s); + return true; + } + catch (Exception) + { + result = 0; + return false; + } + } } } diff --git a/source/Cosmos.System2_Plugs/System/FloatImpl.cs b/source/Cosmos.System2_Plugs/System/FloatImpl.cs new file mode 100644 index 000000000..8ddab07e2 --- /dev/null +++ b/source/Cosmos.System2_Plugs/System/FloatImpl.cs @@ -0,0 +1,118 @@ +using System; +using Cosmos.Common; +using System.Collections.Generic; +using IL2CPU.API.Attribs; + +namespace Cosmos.System_Plugs.System +{ + [Plug(Target = typeof(float))] + public static class FloatImpl + { + public static string ToString(ref float aThis) + { + return StringHelper.GetNumberString(aThis); + } + + public static float Parse(string s) + { + //Format of Float string: [whitespace][sign][integral-digits[,]]integral-digits[.[fractional-digits]][E[sign]exponential-digits][whitespace] + + //Validate input + if (s is null) throw new ArgumentNullException("s can not be null"); + + //Remove leading whitespaces + while (s.Length != 0 && (s[0] == ' ' || s[0] == '\n' || s[0] == '\t')) + { + s = s.Substring(1); + } + + //Check that string is not finished too early + if (s.Length == 0) throw new FormatException(); + + //Check for sign + short sign = 1; + if (s[0] == '-') + { + s = s.Substring(1); + sign = -1; + } + else if (s[0] == '+') s = s.Substring(1); + + //Check that string is not finished too early + if (s.Length == 0) throw new FormatException(); + + //Read in number + + List internalDigits = new List(); + List fractionalDigits = new List(); + + bool foundDecimal = false; + + //Iterate until fully parsed or an E/Whitespace is found + //Assume english standard, so . == decimal seperator and , == thousands seperator + while (s.Length != 0) + { + char active = s[0]; + if (active == 'E' || active == 'e' || active == ' ' || active == '\n' || active == '\t') break; + + s = s.Substring(1); + if (active == '.') foundDecimal = true; + else if (active == ',') continue; + else if (active >= '0' && active <= '9') + { + if (foundDecimal) fractionalDigits.Add(int.Parse(active.ToString())); + else internalDigits.Add(int.Parse(active.ToString())); + } + else + { + throw new FormatException(); + } + } + + //Iterate through rest of string + float multiplier = 0; + while (s.Length != 0) + { + //Check for exponential notation i.e. 8.1E10 = 8.1 * 10^10 + Whitespaces + //E can only be followed by integers + if (s[0] == 'E' || s[0] == 'e') + { + multiplier = float.Parse(s.Substring(1)); + break; + } + else if (s[0] == ' ' || s[0] == '\n' || s[0] == '\t') s = s.Substring(1); + else throw new FormatException(); + } + + //Create float + float parsed = 0; + for (int i = 0; i < internalDigits.Count; i++) + { + parsed += internalDigits[i] * (float)Math.Pow(10, (internalDigits.Count - (i + 1))); + } + for (int i = 0; i < fractionalDigits.Count; i++) + { + parsed += fractionalDigits[i] * (float)Math.Pow(10, -1 * (i + 1)); + } + + parsed *= (float)Math.Pow(10, multiplier); + parsed *= sign; + + return parsed; + } + + public static bool TryParse(string s, out float result) + { + try + { + result = Parse(s); + return true; + } + catch (Exception) + { + result = 0; + return false; + } + } + } +} diff --git a/source/Cosmos.System2_Plugs/System/MathImpl.cs b/source/Cosmos.System2_Plugs/System/MathImpl.cs index f6bf15dfe..33f7842f2 100644 --- a/source/Cosmos.System2_Plugs/System/MathImpl.cs +++ b/source/Cosmos.System2_Plugs/System/MathImpl.cs @@ -11,6 +11,7 @@ namespace Cosmos.System_Plugs.System internal static Debugger mDebugger = new Debugger("System", "Math Plugs"); #region Internal Constants + private const double sq2p1 = 2.414213562373095048802e0F; private const double sq2m1 = .414213562373095048802e0F; private const double pio2 = 1.570796326794896619231e0F; @@ -28,12 +29,14 @@ namespace Cosmos.System_Plugs.System private const double atan_q2 = .16667838148816337184521798e4F; private const double atan_q1 = .207933497444540981287275926e4F; private const double atan_q0 = .89678597403663861962481162e3F; - #endregion + + #endregion Internal Constants public const double PI = 3.1415926535897932384626433832795; public const double E = 2.71828182845904523536; #region Abs + public static double Abs(double value) { if (value < 0) @@ -57,18 +60,22 @@ namespace Cosmos.System_Plugs.System return value; } } - #endregion + + #endregion Abs #region Acos + public static double Acos(double x) { if ((x > 1.0) || (x < -1.0)) throw new ArgumentOutOfRangeException("Domain error"); return (pio2 - Asin(x)); } - #endregion + + #endregion Acos #region Asin + public static double Asin(double x) { if (x > 1.0F) @@ -92,16 +99,20 @@ namespace Cosmos.System_Plugs.System } return (sign * temp); } - #endregion + + #endregion Asin #region Atan + public static double Atan(double x) { return ((x > 0F) ? atans(x) : (-atans(-x))); } - #endregion + + #endregion Atan #region Atan2 + public static double Atan2(double x, double y) { if ((x + y) == x) @@ -126,9 +137,11 @@ namespace Cosmos.System_Plugs.System //return (((x + y) == x) ? (((x == 0F) & (y == 0F)) ? 0F : ((x >= 0F) ? pio2 : (-pio2))) : ((y < 0F) ? ((x >= 0F) ? ((pio2 * 2) - atans((-x) / y)) : (((-pio2) * 2) + atans(x / y))) : ((x > 0F) ? atans(x / y) : -atans((-x) / y)))); } - #endregion + + #endregion Atan2 #region Ceiling + public static double Ceiling(double a) { // should be using assembler for bigger values than int or long max @@ -137,9 +150,11 @@ namespace Cosmos.System_Plugs.System int i = (a - (int)a > 0) ? (int)(a + 1) : (int)a; return i; } - #endregion + + #endregion Ceiling #region Cos + public static double Cos(double x) { // First we need to anchor it to a valid range. @@ -187,55 +202,162 @@ namespace Cosmos.System_Plugs.System return -(c1 + (x2 * (c2 + (x2 * (c3 + (x2 * (c4 + (x2 * (c5 + (x2 * (c6 + (x2 * (c7 + (x2 * (c8 + (x2 * (c9 + (x2 * (c10 + (x2 * c11)))))))))))))))))))); } } - #endregion + + #endregion Cos #region Cosh + public static double Cosh(double x) { if (x < 0.0F) x = -x; return ((x == 0F) ? 1F : ((x <= (ln2 / 2)) ? (1 + (_power((Exp(x) - 1), 2) / (2 * Exp(x)))) : ((x <= 22F) ? ((Exp(x) + (1 / Exp(x))) / 2) : (0.5F * (Exp(x) + Exp(-x)))))); } - #endregion + + #endregion Cosh #region Exp + + /* + * ==================================================== + * Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved. + * + * Permission to use, copy, modify, and distribute this + * software is freely granted, provided that this notice + * is preserved. + * ==================================================== + */ + + // Look at http://www.netlib.org/fdlibm/e_exp.c for more a in deth explanation + private static int HighWord(double x) + { + long value = BitConverter.DoubleToInt64Bits(x); + Byte[] valueBytes = BitConverter.GetBytes(value); + int offset = BitConverter.IsLittleEndian ? 4 : 0; + return BitConverter.ToInt32(valueBytes, offset); + } + + private static int LowWord(double x) //Opposite of high word + { + long value = BitConverter.DoubleToInt64Bits(x); + Byte[] valueBytes = BitConverter.GetBytes(value); + return BitConverter.ToInt32(valueBytes, BitConverter.IsLittleEndian ? 0 : 4); + } + public static double Exp(double x) { - double c; - int n = 1; - double ex = 1F; - double m = 1F; - while (x > 10.000F) - { - m *= 22026.4657948067; - x -= 10F; + double y, hi = 0, lo = 0, c, t; + int k = 0, xsb; + + const double o_threshold = 7.09782712893383973096e+02; + const double u_threshold = -7.45133219101941108420e+02; + const double invln2 = 1.44269504088896338700e+00; + const double twom1000 = 9.33263618503218878990e-302; + const double P1 = 1.66666666666666019037e-01; + const double P2 = -2.77777777770155933842e-03; + const double P3 = 6.61375632143793436117e-05; + const double P4 = -1.65339022054652515390e-06; + const double P5 = 4.13813679705723846039e-08; + const double huge = 1.0e+300; + + int hx = HighWord(x); //Highword of x + + xsb = ((int)hx >> 31) & 1; //Get sign of x + hx &= 0x7fffffff; //Get the abs(x) of the highword + + //Check if non-finite argument + if (hx >= 0x40862E42) + { /* if |x|>=709.78... */ + if (hx >= 0x7ff00000) + { + if (((hx & 0xfffff) | LowWord(x)) != 0) //Assume that __Lo(x) is lower word of x + return x; /* NaN */ + else + return (xsb == 0) ? x : 0.0; /* exp(+-inf)={inf,0} */ + } + if (x > o_threshold) + return double.PositiveInfinity; /* overflow */ + if (x < u_threshold) + return 0; /* underflow */ } - while (x > 01.000F) - { - m *= E; - x -= 1F; + + /* argument reduction */ + if (hx > 0x3fd62e42) + { /* if |x| > 0.5 ln2 */ + if (hx < 0x3FF0A2B2) + { /* and |x| < 1.5 ln2 */ + if (xsb == 0) + { + hi = x - 6.93147180369123816490e-01; + lo = 1.90821492927058770002e-10; + } + else + { + hi = x - -6.93147180369123816490e-01; + lo = -1.90821492927058770002e-10; + } + k = 1 - xsb - xsb; + } + else + { + if (xsb == 0) + k = (int)(invln2 * x + 0.5); + else + k = (int)(invln2 * x + -0.5); + t = k; + hi = x - t * 6.93147180369123816490e-01; + lo = t * 1.90821492927058770002e-10; + } + x = hi - lo; } - while (x > 00.100F) - { - m *= 1.10517091807565; ; - x -= 0.1F; + else if (hx < 0x3e300000) + { /* when |x|<2**-28 */ + if (huge + x > 1) + return 1 + x;/* trigger inexact */ } - while (x > 00.010F) + else + k = 0; + + /* x is now in primary range */ + t = x * x; + c = x - t * (P1 + t * (P2 + t * (P3 + t * (P4 + t * P5)))); + + if (k == 0) + return 1 - ((x * c) / (c - 2.0) - x); + else + y = 1 - ((lo - (x * c) / (2.0 - c)) - hi); + + if (k >= -1021) { - m *= 1.01005016708417; - x -= 0.01F; + //The idea is to add hy to the exponent of y + long _y = BitConverter.DoubleToInt64Bits(y); + + /* add k to y's exponent */ + if (BitConverter.IsLittleEndian) + _y += ((long)k << 52); + else + _y += ((long)k << 20); + y = BitConverter.Int64BitsToDouble(_y); + return y; } - for (int y = 1; y <= 4; y++) + else { - c = _power(x, y); - ex += c / (double)n; - n *= (y + 1); + //The idea is to add hy to the exponent of y + long _y = BitConverter.DoubleToInt64Bits(y); + + if (BitConverter.IsLittleEndian) + _y += ((long)k + 1000 << 52); + else + _y += ((long)k + 1000 << 20); + y = BitConverter.Int64BitsToDouble(_y); + return y * twom1000; } - return ex * m; } - #endregion + + #endregion Exp #region Floor + public static double Floor(double a) { // should be using assembler for bigger values than int or long max @@ -244,16 +366,20 @@ namespace Cosmos.System_Plugs.System int i = (a - (int)a < 0) ? (int)(a - 1) : (int)a; return i; } - #endregion + + #endregion Floor #region Log (base e) + public static double Log(double x) { return Log(x, E); } - #endregion + + #endregion Log (base e) #region Log (base specified) + public static double Log(double x, double newBase) { if (x == 0.0F) @@ -296,65 +422,94 @@ namespace Cosmos.System_Plugs.System return (integer + fractional); } - #endregion + + #endregion Log (base specified) #region Log10 + public static double Log10(double x) { return Log(x, 10F); } - #endregion + + #endregion Log10 #region Pow - public static double Pow(double x, double y) + + public static double Pow(double b, double e) { - if (x <= 0.0F) + if (e == 0) return 1; + if (e == 1) return b; + if (double.IsNaN(b) || double.IsNaN(e)) return double.NaN; + if (double.IsNegativeInfinity(b)) { - double temp = 0F; - long l; - if (x == 0.0F && y <= 0.0F) - throw new ArgumentException(); - - l = (long)Floor(y); - if (l != y) - temp = Exp(y * Log(-x)); - if ((l % 2) == 1) - temp = -temp; - - return (temp); + if (e < 0) + return 0; + if ((long)e % 2 == 0) + return double.PositiveInfinity; + else + return double.NegativeInfinity; } + if (double.IsPositiveInfinity(b)) + { + if (e < 0) + return 0; + else + return double.PositiveInfinity; + } + if (double.IsInfinity(e)) + { + bool t = -1 < b; + bool t1 = 1 > b; + if (t && t1) - return (Exp(y * Log(x))); - //if (y == 0) - //{ - // return 1; - //} - //else if (y == 1) - //{ - // return x; - //} - //else - //{ - // double xResult = x; - - // for (int i = 2; i <= y; i++) - // { - // xResult = xResult * x; - // } - - // return xResult; - //} + if (double.IsPositiveInfinity(e)) + return 0; + else + return double.PositiveInfinity; + else + { + bool v = b < -1; + bool v1 = 1 < b; + if (v || v1) + { + if (double.IsPositiveInfinity(e)) + return double.PositiveInfinity; + else + return 0; + } + else + return double.NaN; + } + } + if (b < 0) + { + if (Math.Abs(e) - Math.Abs((int)e) > (Double.Epsilon * 100)) return double.NaN; + double logedBase = Math.Log(Math.Abs(b)); + double pow = Exp(logedBase * e); + if ((long)e % 2 == 0) return pow; + else return -1 * pow; + } + else + { + double logedBase = Math.Log(b); + return Exp(logedBase * e); + } } - #endregion + + #endregion Pow #region Round + public static double Round(double d) { return ((Math.Floor(d) % 2 == 0) ? Math.Floor(d) : Math.Ceiling(d)); } - #endregion + + #endregion Round #region Sin + public static double Sin(double x) { // First we need to anchor it to a valid range. @@ -365,9 +520,11 @@ namespace Cosmos.System_Plugs.System return Cos((PI / 2.0F) - x); } - #endregion + + #endregion Sin #region Sinh + public static double Sinh(double x) { if (x < 0F) @@ -383,9 +540,11 @@ namespace Cosmos.System_Plugs.System return (Exp(x) / 2); } } - #endregion + + #endregion Sinh #region Sqrt + public static double Sqrt(double x) { long x1; @@ -412,16 +571,18 @@ namespace Cosmos.System_Plugs.System x2 = BitConverter.Int64BitsToDouble(x1); // Use Newton's Method - for(i = 0; i < 5; i++) + for (i = 0; i < 5; i++) { x2 = x2 - (x2 * x2 - x) / (2 * x2); } return x2; } - #endregion + + #endregion Sqrt #region Tan + public static double Tan(double x) { // First we need to anchor it to a valid range. @@ -436,30 +597,36 @@ namespace Cosmos.System_Plugs.System case 0: x = x * (4 / PI); break; + case 1: x = ((PI / 2) - x) * (4 / PI); break; + case 2: x = (x - (PI / 2)) * (4 / PI); break; + case 3: x = (PI - x) * (4 / PI); break; + case 4: x = (x - PI) * (4 / PI); break; + case 5: x = ((3.5 * PI) - x) * (4 / PI); break; + case 6: x = (x - (3.5 * PI)) * (4 / PI); break; + case 7: x = ((2 * PI) - x) * (4 / PI); break; } - const double c1 = 4130240.588996024013440146267; const double c2 = -349781.8562517381616631012487; const double c3 = 6170.317758142494245331944348; @@ -487,21 +654,26 @@ namespace Cosmos.System_Plugs.System return -((x * (c1 + (x2 * (c2 + (x2 * (c3 + (x2 * (c4 + (x2 * c5))))))))) / (c6 + (x2 * (c7 + (x2 * (c8 + (x2 * (c9 + x2)))))))); } } - #endregion + + #endregion Tan #region Tanh + public static double Tanh(double x) { return (expm1(2F * x) / (expm1(2F * x) + 2F)); } - #endregion + + #endregion Tanh #region Truncate + public static double Truncate(double x) { return ((x == 0) ? 0F : ((x > 0F) ? Floor(x) : Ceiling(x))); } - #endregion + + #endregion Truncate //#region Factorial (only used in Sin(), not plug ) //public static int Factorial(int n) @@ -516,14 +688,17 @@ namespace Cosmos.System_Plugs.System #region Internaly used functions #region expm1 + private static double expm1(double x) { double u = Exp(x); return ((u == 1.0F) ? x : ((u - 1.0F == -1.0F) ? -1.0F : ((u - 1.0F) * x / Log(u)))); } - #endregion + + #endregion expm1 #region _power + private static double _power(double x, int c) { if (c == 0) @@ -545,9 +720,11 @@ namespace Cosmos.System_Plugs.System return ret; } - #endregion + + #endregion _power #region atans + private static double atans(double x) { if (x < sq2m1) @@ -563,9 +740,11 @@ namespace Cosmos.System_Plugs.System return (pio4 + atanx((x - 1.0F) / (x + 1.0F))); } } - #endregion + + #endregion atans #region atanx + private static double atanx(double x) { double argsq, value; @@ -592,8 +771,9 @@ namespace Cosmos.System_Plugs.System // / // (((((ArgSquared + atan_q4) * ArgSquared + atan_q3) * ArgSquared + atan_q2) * ArgSquared + atan_q1) * ArgSquared + atan_q0) * x); } - #endregion - #endregion + #endregion atanx + + #endregion Internaly used functions } }