diff --git a/src/Numerics.Tests/IntegrationTests/IntegrationTest.cs b/src/Numerics.Tests/IntegrationTests/IntegrationTest.cs index dd988098..bc2e967b 100644 --- a/src/Numerics.Tests/IntegrationTests/IntegrationTest.cs +++ b/src/Numerics.Tests/IntegrationTests/IntegrationTest.cs @@ -181,16 +181,13 @@ namespace MathNet.Numerics.UnitTests.IntegrationTests TargetAreaC, Integrate.DoubleExponential(TargetFunctionC, StartC, StopC, 1e-10), 1e-10, - "DoubleExponential, Target 1e-10"); - - // integrate_(0)^(1) log(x) dx = -1 - // Note that DoubleExponential returns -oo. + "DoubleExponential"); Assert.AreEqual( TargetAreaD, - Integrate.OnClosedInterval(TargetFunctionD, StartD, StopD), + Integrate.DoubleExponential(TargetFunctionD, StartD, StopD), 1e-10, - "Interval"); + "DoubleExponential"); Assert.AreEqual( TargetAreaD, @@ -474,19 +471,19 @@ namespace MathNet.Numerics.UnitTests.IntegrationTests expected, factor * Integrate.DoubleExponential((x) => 1 / (1 + x * x), a, b), 1e-10, - "DET Integral of sin(pi*x)/(pi*x) from -oo to oo"); + "DET Integral of sin(pi*x)/(pi*x) from {0} to {1}", a, b); Assert.AreEqual( expected, factor * Integrate.GaussKronrod((x) => 1 / (1 + x * x), a, b), 1e-10, - "GK Integral of sin(pi*x)/(pi*x) from -oo to oo"); + "GK Integral of sin(pi*x)/(pi*x) from {0} to {1}", a, b); Assert.AreEqual( expected, factor * Integrate.GaussLegendre((x) => 1 / (1 + x * x), a, b, order: 128), 1e-10, - "GL Integral of sin(pi*x)/(pi*x) from -oo to oo"); + "GL Integral of sin(pi*x)/(pi*x) from {0} to {1}", a, b); } // integral_(-oo)^(oo) 1/(1 + j x^2) dx = -(-1)^(3/4) ¥ð diff --git a/src/Numerics/Integrate.cs b/src/Numerics/Integrate.cs index 9c7e0370..5761f6c9 100644 --- a/src/Numerics/Integrate.cs +++ b/src/Numerics/Integrate.cs @@ -102,7 +102,7 @@ namespace MathNet.Numerics // (-oo, oo) => [-1, 1] // - // integral_{-oo}^{oo} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(-oo)^(oo) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = t / (1 - t^2) // g'(t) = (1 + t^2) / (1 - t^2)^2 if (double.IsInfinity(intervalBegin) && double.IsInfinity(intervalEnd)) @@ -115,8 +115,8 @@ namespace MathNet.Numerics } // [a, oo) => [0, 1] // - // integral_{a}^{oo} f(x) dx = integral_{0}^{oo} f(a + t^2) 2 t dt - // = integral_{0}^{1} f(a + g(s)^2) 2 g(s) g'(s) ds + // integral_(a)^(oo) f(x) dx = integral_(0)^(oo) f(a + t^2) 2 t dt + // = integral_(0)^(1) f(a + g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 - s) // g'(s) = 1 / (1 - s)^2 else if (double.IsInfinity(intervalEnd)) @@ -129,8 +129,8 @@ namespace MathNet.Numerics } // (-oo, b] => [-1, 0] // - // integral_{-oo}^{b} f(x) dx = -integral_{-oo}^{0} f(b - t^2) 2 t dt - // = -integral_{-1}^{0} f(b - g(s)^2) 2 g(s) g'(s) ds + // integral_(-oo)^(b) f(x) dx = -integral_(-oo)^(0) f(b - t^2) 2 t dt + // = -integral_(-1)^(0) f(b - g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 + s) // g'(s) = 1 / (1 + s)^2 else if (double.IsInfinity(intervalBegin)) @@ -141,18 +141,9 @@ namespace MathNet.Numerics }; return DoubleExponentialTransformation.Integrate(u, -1, 0, targetAbsoluteError); } - // [a, b] => [-1, 1] - // - // integral_{a}^{b} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt - // g(t) = (b - a) * t * (3 - t^2) / 4 + (b + a) / 2 - // g'(t) = 3 / 4 * (b - a) * (1 - t^2) else { - Func u = (t) => - { - return f((intervalEnd - intervalBegin) / 4 * t * (3 - t * t) + (intervalEnd + intervalBegin) / 2) * 3 * (intervalEnd - intervalBegin) / 4 * (1 - t * t); - }; - return DoubleExponentialTransformation.Integrate(u, -1, 1, targetAbsoluteError); + return DoubleExponentialTransformation.Integrate(f, intervalBegin, intervalEnd, targetAbsoluteError); } } @@ -178,7 +169,7 @@ namespace MathNet.Numerics // (-oo, oo) => [-1, 1] // - // integral_{-oo}^{oo} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(-oo)^(oo) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = t / (1 - t^2) // g'(t) = (1 + t^2) / (1 - t^2)^2 if (double.IsInfinity(intervalBegin) && double.IsInfinity(intervalEnd)) @@ -191,8 +182,8 @@ namespace MathNet.Numerics } // [a, oo) => [0, 1] // - // integral_{a}^{oo} f(x) dx = integral_{0}^{oo} f(a + t^2) 2 t dt - // = integral_{0}^{1} f(a + g(s)^2) 2 g(s) g'(s) ds + // integral_(a)^(oo) f(x) dx = integral_(0)^(oo) f(a + t^2) 2 t dt + // = integral_(0)^(1) f(a + g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 - s) // g'(s) = 1 / (1 - s)^2 else if (double.IsInfinity(intervalEnd)) @@ -205,8 +196,8 @@ namespace MathNet.Numerics } // (-oo, b] => [-1, 0] // - // integral_{-oo}^{b} f(x) dx = -integral_{-oo}^{0} f(b - t^2) 2 t dt - // = -integral_{-1}^{0} f(b - g(s)^2) 2 g(s) g'(s) ds + // integral_(-oo)^(b) f(x) dx = -integral_(-oo)^(0) f(b - t^2) 2 t dt + // = -integral_(-1)^(0) f(b - g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 + s) // g'(s) = 1 / (1 + s)^2 else if (double.IsInfinity(intervalBegin)) @@ -219,7 +210,7 @@ namespace MathNet.Numerics } // [a, b] => [-1, 1] // - // integral_{a}^{b} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(a)^(b) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = (b - a) * t * (3 - t^2) / 4 + (b + a) / 2 // g'(t) = 3 / 4 * (b - a) * (1 - t^2) else @@ -292,7 +283,7 @@ namespace MathNet.Numerics // (-oo, oo) => [-1, 1] // - // integral_{-oo}^{oo} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(-oo)^(oo) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = t / (1 - t^2) // g'(t) = (1 + t^2) / (1 - t^2)^2 if (double.IsInfinity(intervalBegin) && double.IsInfinity(intervalEnd)) @@ -305,8 +296,8 @@ namespace MathNet.Numerics } // [a, oo) => [0, 1] // - // integral_{a}^{oo} f(x) dx = integral_{0}^{oo} f(a + t^2) 2 t dt - // = integral_{0}^{1} f(a + g(s)^2) 2 g(s) g'(s) ds + // integral_(a)^(oo) f(x) dx = integral_(0)^(oo) f(a + t^2) 2 t dt + // = integral_(0)^(1) f(a + g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 - s) // g'(s) = 1 / (1 - s)^2 else if (double.IsInfinity(intervalEnd)) @@ -319,8 +310,8 @@ namespace MathNet.Numerics } // (-oo, b] => [-1, 0] // - // integral_{-oo}^{b} f(x) dx = -integral_{-oo}^{0} f(b - t^2) 2 t dt - // = -integral_{-1}^{0} f(b - g(s)^2) 2 g(s) g'(s) ds + // integral_(-oo)^(b) f(x) dx = -integral_(-oo)^(0) f(b - t^2) 2 t dt + // = -integral_(-1)^(0) f(b - g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 + s) // g'(s) = 1 / (1 + s)^2 else if (double.IsInfinity(intervalBegin)) @@ -331,18 +322,9 @@ namespace MathNet.Numerics }; return DoubleExponentialTransformation.ContourIntegrate(u, -1, 0, targetAbsoluteError); } - // [a, b] => [-1, 1] - // - // integral_{a}^{b} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt - // g(t) = (b - a) * t * (3 - t^2) / 4 + (b + a) / 2 - // g'(t) = 3 / 4 * (b - a) * (1 - t^2) else { - Func u = (t) => - { - return f((intervalEnd - intervalBegin) / 4 * t * (3 - t * t) + (intervalEnd + intervalBegin) / 2) * 3 * (intervalEnd - intervalBegin) / 4 * (1 - t * t); - }; - return DoubleExponentialTransformation.ContourIntegrate(u, -1, 1, targetAbsoluteError); + return DoubleExponentialTransformation.ContourIntegrate(f, intervalBegin, intervalEnd, targetAbsoluteError); } } @@ -368,7 +350,7 @@ namespace MathNet.Numerics // (-oo, oo) => [-1, 1] // - // integral_{-oo}^{oo} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(-oo)^(oo) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = t / (1 - t^2) // g'(t) = (1 + t^2) / (1 - t^2)^2 if (double.IsInfinity(intervalBegin) && double.IsInfinity(intervalEnd)) @@ -381,8 +363,8 @@ namespace MathNet.Numerics } // [a, oo) => [0, 1] // - // integral_{a}^{oo} f(x) dx = integral_{0}^{oo} f(a + t^2) 2 t dt - // = integral_{0}^{1} f(a + g(s)^2) 2 g(s) g'(s) ds + // integral_(a)^(oo) f(x) dx = integral_(0)^(oo) f(a + t^2) 2 t dt + // = integral_(0)^(1) f(a + g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 - s) // g'(s) = 1 / (1 - s)^2 else if (double.IsInfinity(intervalEnd)) @@ -395,8 +377,8 @@ namespace MathNet.Numerics } // (-oo, b] => [-1, 0] // - // integral_{-oo}^{b} f(x) dx = -integral_{-oo}^{0} f(b - t^2) 2 t dt - // = -integral_{-1}^{0} f(b - g(s)^2) 2 g(s) g'(s) ds + // integral_(-oo)^(b) f(x) dx = -integral_(-oo)^(0) f(b - t^2) 2 t dt + // = -integral_(-1)^(0) f(b - g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 + s) // g'(s) = 1 / (1 + s)^2 else if (double.IsInfinity(intervalBegin)) @@ -409,7 +391,7 @@ namespace MathNet.Numerics } // [a, b] => [-1, 1] // - // integral_{a}^{b} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(a)^(b) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = (b - a) * t * (3 - t^2) / 4 + (b + a) / 2 // g'(t) = 3 / 4 * (b - a) * (1 - t^2) else diff --git a/src/Numerics/Integration/GaussKronrodRule.cs b/src/Numerics/Integration/GaussKronrodRule.cs index 82e1bf36..3d5224be 100644 --- a/src/Numerics/Integration/GaussKronrodRule.cs +++ b/src/Numerics/Integration/GaussKronrodRule.cs @@ -240,7 +240,7 @@ namespace MathNet.Numerics.Integration // (-oo, oo) => [-1, 1] // - // integral_{-oo}^{oo} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(-oo)^(oo) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = t / (1 - t^2) // g'(t) = (1 + t^2) / (1 - t^2)^2 if ((intervalBegin < double.MinValue) && (intervalEnd > double.MaxValue)) @@ -253,8 +253,8 @@ namespace MathNet.Numerics.Integration } // [a, oo) => [0, 1] // - // integral_{a}^{oo} f(x) dx = integral_{0}^{oo} f(a + t^2) 2 t dt - // = integral_{0}^{1} f(a + g(s)^2) 2 g(s) g'(s) ds + // integral_(a)^(oo) f(x) dx = integral_(0)^(oo) f(a + t^2) 2 t dt + // = integral_(0)^(1) f(a + g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 - s) // g'(s) = 1 / (1 - s)^2 else if (intervalEnd > double.MaxValue) @@ -267,8 +267,8 @@ namespace MathNet.Numerics.Integration } // (-oo, b] => [-1, 0] // - // integral_{-oo}^{b} f(x) dx = -integral_{-oo}^{0} f(b - t^2) 2 t dt - // = -integral_{-1}^{0} f(b - g(s)^2) 2 g(s) g'(s) ds + // integral_(-oo)^(b) f(x) dx = -integral_(-oo)^(0) f(b - t^2) 2 t dt + // = -integral_(-1)^(0) f(b - g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 + s) // g'(s) = 1 / (1 + s)^2 else if (intervalBegin < double.MinValue) @@ -281,7 +281,7 @@ namespace MathNet.Numerics.Integration } // [a, b] => [-1, 1] // - // integral_{a}^{b} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(a)^(b) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = (b - a) * t * (3 - t^2) / 4 + (b + a) / 2 // g'(t) = 3 / 4 * (b - a) * (1 - t^2) else @@ -326,7 +326,7 @@ namespace MathNet.Numerics.Integration // (-oo, oo) => [-1, 1] // - // integral_{-oo}^{oo} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(-oo)^(oo) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = t / (1 - t^2) // g'(t) = (1 + t^2) / (1 - t^2)^2 if ((intervalBegin < double.MinValue) && (intervalEnd > double.MaxValue)) @@ -339,8 +339,8 @@ namespace MathNet.Numerics.Integration } // [a, oo) => [0, 1] // - // integral_{a}^{oo} f(x) dx = integral_{0}^{oo} f(a + t^2) 2 t dt - // = integral_{0}^{1} f(a + g(s)^2) 2 g(s) g'(s) ds + // integral_(a)^(oo) f(x) dx = integral_(0)^(oo) f(a + t^2) 2 t dt + // = integral_(0)^(1) f(a + g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 - s) // g'(s) = 1 / (1 - s)^2 else if (intervalEnd > double.MaxValue) @@ -353,8 +353,8 @@ namespace MathNet.Numerics.Integration } // (-oo, b] => [-1, 0] // - // integral_{-oo}^{b} f(x) dx = -integral_{-oo}^{0} f(b - t^2) 2 t dt - // = -integral_{-1}^{0} f(b - g(s)^2) 2 g(s) g'(s) ds + // integral_(-oo)^(b) f(x) dx = -integral_(-oo)^(0) f(b - t^2) 2 t dt + // = -integral_(-1)^(0) f(b - g(s)^2) 2 g(s) g'(s) ds // g(s) = s / (1 + s) // g'(s) = 1 / (1 + s)^2 else if (intervalBegin < double.MinValue) @@ -367,7 +367,7 @@ namespace MathNet.Numerics.Integration } // [a, b] => [-1, 1] // - // integral_{a}^{b} f(x) dx = integral_{-1}^{1} f(g(t)) g'(t) dt + // integral_(a)^(b) f(x) dx = integral_(-1)^(1) f(g(t)) g'(t) dt // g(t) = (b - a) * t * (3 - t^2) / 4 + (b + a) / 2 // g'(t) = 3 / 4 * (b - a) * (1 - t^2) else