From 06730b30d64b6012ddebf0c6c5958cd0b3aa0830 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?G=C3=BCnther=20Foidl?= Date: Thu, 19 Apr 2018 18:24:18 +0200 Subject: [PATCH 1/2] Added tests for Sample with size 101 To test uneven size, especially for loop-unrolling. --- tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Average.cs | 1 + tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Delta.cs | 1 + tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs | 1 + .../Statistics/SampleTests/SkewnessAndKurtosis.cs | 1 + .../Statistics/SampleTests/VarianceCore.cs | 1 + .../Statistics/SampleTests/ZTransformation.cs | 1 + 6 files changed, 6 insertions(+) diff --git a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Average.cs b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Average.cs index 7e646f4..aaa51bd 100644 --- a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Average.cs +++ b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Average.cs @@ -7,6 +7,7 @@ namespace gfoidl.Stochastics.Tests.Statistics.SampleTests { [TestFixture(10)] [TestFixture(100)] + [TestFixture(101)] [TestFixture(1_000)] [TestFixture(10_000)] [TestFixture(100_000)] diff --git a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Delta.cs b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Delta.cs index d182ec9..bdd47fe 100644 --- a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Delta.cs +++ b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/Delta.cs @@ -6,6 +6,7 @@ namespace gfoidl.Stochastics.Tests.Statistics.SampleTests { [TestFixture(10)] [TestFixture(100)] + [TestFixture(101)] [TestFixture(1_000)] [TestFixture(10_000)] [TestFixture(100_000)] diff --git a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs index fa56395..96f4434 100644 --- a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs +++ b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs @@ -7,6 +7,7 @@ namespace gfoidl.Stochastics.Tests.Statistics.SampleTests { [TestFixture(10)] [TestFixture(100)] + [TestFixture(101)] [TestFixture(1_000)] [TestFixture(10_000)] [TestFixture(100_000)] diff --git a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/SkewnessAndKurtosis.cs b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/SkewnessAndKurtosis.cs index 852a596..b53ec8b 100644 --- a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/SkewnessAndKurtosis.cs +++ b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/SkewnessAndKurtosis.cs @@ -6,6 +6,7 @@ namespace gfoidl.Stochastics.Tests.Statistics.SampleTests { [TestFixture(10)] [TestFixture(100)] + [TestFixture(101)] [TestFixture(1_000)] [TestFixture(10_000)] [TestFixture(100_000)] diff --git a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/VarianceCore.cs b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/VarianceCore.cs index 97637b6..6919a6f 100644 --- a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/VarianceCore.cs +++ b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/VarianceCore.cs @@ -6,6 +6,7 @@ namespace gfoidl.Stochastics.Tests.Statistics.SampleTests { [TestFixture(10)] [TestFixture(100)] + [TestFixture(101)] [TestFixture(1_000)] [TestFixture(10_000)] [TestFixture(100_000)] diff --git a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/ZTransformation.cs b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/ZTransformation.cs index 465f609..a1932a2 100644 --- a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/ZTransformation.cs +++ b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/ZTransformation.cs @@ -7,6 +7,7 @@ namespace gfoidl.Stochastics.Tests.Statistics.SampleTests { [TestFixture(10)] [TestFixture(100)] + [TestFixture(101)] [TestFixture(1_000)] [TestFixture(10_000)] [TestFixture(100_000)] From 6a2f9d0db1809c45ee7a3b43f7fdee329c8d39cc Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?G=C3=BCnther=20Foidl?= Date: Thu, 19 Apr 2018 18:12:42 +0200 Subject: [PATCH 2/2] Loop limit determined by closest multiple (bit hack) instead of naive offset --- .../Statistics/Sample.AutoCorrelation.cs | 10 +-- .../Statistics/Sample.AverageVarianceCore.cs | 70 +++++++++------ .../Statistics/Sample.Delta.cs | 70 +++++++++------ .../Statistics/Sample.MinMax.cs | 66 ++++++++------ .../Statistics/Sample.SkewnessAndKurtosis.cs | 76 ++++++++++------- .../Statistics/Sample.ZTransformation.cs | 85 +++++++++++-------- .../gfoidl.Stochastics.csproj | 6 +- .../Statistics/SampleTests/MinMax.cs | 22 +++-- 8 files changed, 251 insertions(+), 154 deletions(-) diff --git a/source/gfoidl.Stochastics/Statistics/Sample.AutoCorrelation.cs b/source/gfoidl.Stochastics/Statistics/Sample.AutoCorrelation.cs index bf5f75b..a197fce 100644 --- a/source/gfoidl.Stochastics/Statistics/Sample.AutoCorrelation.cs +++ b/source/gfoidl.Stochastics/Statistics/Sample.AutoCorrelation.cs @@ -38,13 +38,13 @@ internal IEnumerable AutoCorrelationSimd() { var kVec = new Vector(arr, k); var kmVec = new Vector(arr, k - m); - r_xx += Vector.Dot(kVec, kmVec); - k += Vector.Count; + r_xx += Vector.Dot(kVec, kmVec); + k += Vector.Count; kVec = new Vector(arr, k); kmVec = new Vector(arr, k - m); r_xx += Vector.Dot(kVec, kmVec); - k += Vector.Count; + k += Vector.Count; } for (; k < arr.Length; ++k) @@ -84,7 +84,7 @@ internal double[] AutoCorrelationToArrayParallelSimd() return corr; } //--------------------------------------------------------------------- - private unsafe void AutoCorrelationToArrayImpl(double[] corr, (int Start, int End) range) + private unsafe void AutoCorrelationToArrayImpl(double[] corr, Range range) { int n = _values.Length; @@ -157,7 +157,7 @@ void Core(double* a_k, double* a_km, int offset, ref double r_xx) { Vector kVec = VectorHelper.GetVector(a_k + offset); Vector kmVec = VectorHelper.GetVector(a_km + offset); - r_xx += Vector.Dot(kVec, kmVec); + r_xx += Vector.Dot(kVec, kmVec); } } #endregion diff --git a/source/gfoidl.Stochastics/Statistics/Sample.AverageVarianceCore.cs b/source/gfoidl.Stochastics/Statistics/Sample.AverageVarianceCore.cs index 35f1b5d..5fcb8b4 100644 --- a/source/gfoidl.Stochastics/Statistics/Sample.AverageVarianceCore.cs +++ b/source/gfoidl.Stochastics/Statistics/Sample.AverageVarianceCore.cs @@ -2,6 +2,10 @@ using System.Numerics; using System.Threading.Tasks; +#if DEBUG_ASSERT +using System.Diagnostics; +#endif + namespace gfoidl.Stochastics.Statistics { partial class Sample @@ -55,71 +59,83 @@ private unsafe void CalculateAverageAndVarianceCoreImpl(int i, int n, out double fixed (double* pArray = _values) { double* arr = pArray + i; + n -= i; + i = 0; + double* end = arr + n; - if (Vector.IsHardwareAccelerated && (n - i) >= Vector.Count) + if (Vector.IsHardwareAccelerated && n >= Vector.Count) { var avgVec = Vector.Zero; - for (; i < n - 8 * Vector.Count; i += 8 * Vector.Count) + // https://github.com/gfoidl/Stochastics/issues/46 + int m = n & ~(8 * Vector.Count - 1); + for (; i < m; i += 8 * Vector.Count) { - Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 1 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 2 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 3 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 4 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 5 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 6 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 7 * Vector.Count, ref avgVec, ref tmpVariance); + Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 1 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 2 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 3 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 4 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 5 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 6 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 7 * Vector.Count, ref avgVec, ref tmpVariance, end); arr += 8 * Vector.Count; } - if (i < n - 4 * Vector.Count) + m = n & ~(4 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 1 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 2 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 3 * Vector.Count, ref avgVec, ref tmpVariance); + Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 1 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 2 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 3 * Vector.Count, ref avgVec, ref tmpVariance, end); arr += 4 * Vector.Count; i += 4 * Vector.Count; } - if (i < n - 2 * Vector.Count) + m = n & ~(2 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance); - Core(arr, 1 * Vector.Count, ref avgVec, ref tmpVariance); + Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance, end); + Core(arr, 1 * Vector.Count, ref avgVec, ref tmpVariance, end); arr += 2 * Vector.Count; i += 2 * Vector.Count; } - if (i < n - Vector.Count) + m = n & ~(1 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance); + Core(arr, 0 * Vector.Count, ref avgVec, ref tmpVariance, end); - i += Vector.Count; + arr += 1 * Vector.Count; } // Reduction -- https://github.com/gfoidl/Stochastics/issues/43 tmpAvg += avgVec.ReduceSum(); } - for (; i < n; ++i) + while (arr < end) { - tmpAvg += pArray[i]; - tmpVariance += pArray[i] * pArray[i]; + tmpAvg += *arr; + tmpVariance += *arr * *arr; + arr++; } avg = tmpAvg; variance = tmpVariance; } //----------------------------------------------------------------- - void Core(double* arr, int offset, ref Vector avgVec, ref double var) + void Core(double* arr, int offset, ref Vector avgVec, ref double var, double* end) { +#if DEBUG_ASSERT + Debug.Assert(arr + offset < end); +#endif Vector vec = VectorHelper.GetVector(arr + offset); - avgVec += vec; - var += Vector.Dot(vec, vec); + avgVec += vec; + var += Vector.Dot(vec, vec); } } } diff --git a/source/gfoidl.Stochastics/Statistics/Sample.Delta.cs b/source/gfoidl.Stochastics/Statistics/Sample.Delta.cs index 25025d2..3942a91 100644 --- a/source/gfoidl.Stochastics/Statistics/Sample.Delta.cs +++ b/source/gfoidl.Stochastics/Statistics/Sample.Delta.cs @@ -3,6 +3,10 @@ using System.Numerics; using System.Threading.Tasks; +#if DEBUG_ASSERT +using System.Diagnostics; +#endif + namespace gfoidl.Stochastics.Statistics { partial class Sample @@ -45,67 +49,81 @@ private unsafe double CalculateDeltaImpl(int i, int n) fixed (double* pArray = _values) { double* arr = pArray + i; + n -= i; + i = 0; + double* end = arr + n; - if (Vector.IsHardwareAccelerated && (n - i) >= Vector.Count) + if (Vector.IsHardwareAccelerated && n >= Vector.Count) { var avgVec = new Vector(avg); - var deltaVec = new Vector(0); + var deltaVec = Vector.Zero; - for (; i < n - 8 * Vector.Count; i += 8 * Vector.Count) + // https://github.com/gfoidl/Stochastics/issues/46 + int m = n & ~(8 * Vector.Count - 1); + for (; i < m; i += 8 * Vector.Count) { - Core(arr, 0 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 1 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 2 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 3 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 4 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 5 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 6 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 7 * Vector.Count, avgVec, ref deltaVec); + Core(arr, 0 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 1 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 2 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 3 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 4 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 5 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 6 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 7 * Vector.Count, avgVec, ref deltaVec, end); arr += 8 * Vector.Count; } - if (i < n - 4 * Vector.Count) + m = n & ~(4 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 1 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 2 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 3 * Vector.Count, avgVec, ref deltaVec); + Core(arr, 0 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 1 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 2 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 3 * Vector.Count, avgVec, ref deltaVec, end); arr += 4 * Vector.Count; i += 4 * Vector.Count; } - if (i < n - 2 * Vector.Count) + m = n & ~(2 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, avgVec, ref deltaVec); - Core(arr, 1 * Vector.Count, avgVec, ref deltaVec); + Core(arr, 0 * Vector.Count, avgVec, ref deltaVec, end); + Core(arr, 1 * Vector.Count, avgVec, ref deltaVec, end); arr += 2 * Vector.Count; i += 2 * Vector.Count; } - if (i < n - Vector.Count) + m = n & ~(1 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, avgVec, ref deltaVec); + Core(arr, 0 * Vector.Count, avgVec, ref deltaVec, end); - i += Vector.Count; + arr += 1 * Vector.Count; } // Reduction -- https://github.com/gfoidl/Stochastics/issues/43 delta += deltaVec.ReduceSum(); } - for (; i < n; ++i) - delta += Math.Abs(pArray[i] - avg); + while (arr < end) + { + delta += Math.Abs(*arr - avg); + arr++; + } } return delta; //----------------------------------------------------------------- - void Core(double* arr, int offset, Vector avgVec, ref Vector deltaVec) + void Core(double* arr, int offset, Vector avgVec, ref Vector deltaVec, double* end) { +#if DEBUG_ASSERT + Debug.Assert(arr + offset < end); +#endif Vector vec = VectorHelper.GetVector(arr + offset); - deltaVec += Vector.Abs(vec - avgVec); + deltaVec += Vector.Abs(vec - avgVec); } } } diff --git a/source/gfoidl.Stochastics/Statistics/Sample.MinMax.cs b/source/gfoidl.Stochastics/Statistics/Sample.MinMax.cs index 716acff..7ef524e 100644 --- a/source/gfoidl.Stochastics/Statistics/Sample.MinMax.cs +++ b/source/gfoidl.Stochastics/Statistics/Sample.MinMax.cs @@ -2,6 +2,10 @@ using System.Numerics; using System.Threading.Tasks; +#if DEBUG_ASSERT +using System.Diagnostics; +#endif + namespace gfoidl.Stochastics.Statistics { partial class Sample @@ -54,69 +58,81 @@ private unsafe void GetMinMaxImpl(int i, int n, out double min, out double max) fixed (double* pArray = _values) { double* arr = pArray + i; + n -= i; + i = 0; + double* end = arr + n; - if (Vector.IsHardwareAccelerated && (n - i) >= Vector.Count) + if (Vector.IsHardwareAccelerated && n >= Vector.Count) { var minVec = new Vector(tmpMin); var maxVec = new Vector(tmpMax); - for (; i < n - 8 * Vector.Count; i += 8 * Vector.Count) + // https://github.com/gfoidl/Stochastics/issues/46 + int m = n & ~(8 * Vector.Count - 1); + for (; i < m; i += 8 * Vector.Count) { - Core(arr, 0 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 1 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 2 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 3 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 4 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 5 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 6 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 7 * Vector.Count, ref minVec, ref maxVec); + Core(arr, 0 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 1 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 2 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 3 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 4 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 5 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 6 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 7 * Vector.Count, ref minVec, ref maxVec, end); arr += 8 * Vector.Count; } - if (i < n - 4 * Vector.Count) + m = n & ~(4 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 1 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 2 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 3 * Vector.Count, ref minVec, ref maxVec); + Core(arr, 0 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 1 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 2 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 3 * Vector.Count, ref minVec, ref maxVec, end); arr += 4 * Vector.Count; i += 4 * Vector.Count; } - if (i < n - 2 * Vector.Count) + m = n & ~(2 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, ref minVec, ref maxVec); - Core(arr, 1 * Vector.Count, ref minVec, ref maxVec); + Core(arr, 0 * Vector.Count, ref minVec, ref maxVec, end); + Core(arr, 1 * Vector.Count, ref minVec, ref maxVec, end); arr += 2 * Vector.Count; i += 2 * Vector.Count; } - if (i < n - Vector.Count) + m = n & ~(1 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, ref minVec, ref maxVec); + Core(arr, 0 * Vector.Count, ref minVec, ref maxVec, end); - i += Vector.Count; + arr += 1 * Vector.Count; } // Reduction VectorHelper.ReduceMinMax(minVec, maxVec, ref tmpMin, ref tmpMax); } - for (; i < n; ++i) + while (arr < end) { - if (pArray[i] < tmpMin) tmpMin = pArray[i]; - if (pArray[i] > tmpMax) tmpMax = pArray[i]; + if (*arr < tmpMin) tmpMin = *arr; + if (*arr > tmpMax) tmpMax = *arr; + arr++; } min = tmpMin; max = tmpMax; } //----------------------------------------------------------------- - void Core(double* arr, int offset, ref Vector minVec, ref Vector maxVec) + void Core(double* arr, int offset, ref Vector minVec, ref Vector maxVec, double* end) { +#if DEBUG_ASSERT + Debug.Assert(arr + offset < end); +#endif Vector vec = VectorHelper.GetVector(arr + offset); minVec = Vector.Min(minVec, vec); maxVec = Vector.Max(maxVec, vec); diff --git a/source/gfoidl.Stochastics/Statistics/Sample.SkewnessAndKurtosis.cs b/source/gfoidl.Stochastics/Statistics/Sample.SkewnessAndKurtosis.cs index bb71b96..301fb8d 100644 --- a/source/gfoidl.Stochastics/Statistics/Sample.SkewnessAndKurtosis.cs +++ b/source/gfoidl.Stochastics/Statistics/Sample.SkewnessAndKurtosis.cs @@ -2,6 +2,10 @@ using System.Numerics; using System.Threading.Tasks; +#if DEBUG_ASSERT +using System.Diagnostics; +#endif + namespace gfoidl.Stochastics.Statistics { partial class Sample @@ -58,52 +62,60 @@ private unsafe void CalculateSkewnessAndKurtosisImpl(int i, int n, out double sk fixed (double* pArray = _values) { double* arr = pArray + i; + n -= i; + i = 0; + double* end = arr + n; - if (Vector.IsHardwareAccelerated && (n - i) >= Vector.Count) + if (Vector.IsHardwareAccelerated && n >= Vector.Count) { var avgVec = new Vector(avg); - var skewVec = new Vector(0); - var kurtVec = new Vector(0); + var skewVec = Vector.Zero; + var kurtVec = Vector.Zero; - for (; i < n - 8 * Vector.Count; i += 8 * Vector.Count) + // https://github.com/gfoidl/Stochastics/issues/46 + int m = n & ~(8 * Vector.Count - 1); + for (; i < m; i += 8 * Vector.Count) { - Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 1 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 2 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 3 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 4 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 5 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 6 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 7 * Vector.Count, avgVec, ref skewVec, ref kurtVec); + Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 1 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 2 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 3 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 4 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 5 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 6 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 7 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); arr += 8 * Vector.Count; } - if (i < n - 4 * Vector.Count) + m = n & ~(4 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 1 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 2 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 3 * Vector.Count, avgVec, ref skewVec, ref kurtVec); + Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 1 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 2 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 3 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); arr += 4 * Vector.Count; i += 4 * Vector.Count; } - if (i < n - 2 * Vector.Count) + m = n & ~(2 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec); - Core(arr, 1 * Vector.Count, avgVec, ref skewVec, ref kurtVec); + Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); + Core(arr, 1 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); arr += 2 * Vector.Count; i += 2 * Vector.Count; } - if (i < n - Vector.Count) + m = n & ~(1 * Vector.Count - 1); + if (i < m) { - Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec); + Core(arr, 0 * Vector.Count, avgVec, ref skewVec, ref kurtVec, end); - i += Vector.Count; + arr += 1 * Vector.Count; } // Reduction -- https://github.com/gfoidl/Stochastics/issues/43 @@ -111,25 +123,29 @@ private unsafe void CalculateSkewnessAndKurtosisImpl(int i, int n, out double sk tmpKurtosis += kurtVec.ReduceSum(); } - for (; i < n; ++i) + while (arr < end) { - double t = pArray[i] - avg; - double t1 = t * t * t; + double t = *arr - avg; + double t1 = t * t * t; tmpSkewness += t1; tmpKurtosis += t1 * t; + arr++; } skewness = tmpSkewness; kurtosis = tmpKurtosis; } //----------------------------------------------------------------- - void Core(double* arr, int offset, Vector avgVec, ref Vector skewVec, ref Vector kurtVec) + void Core(double* arr, int offset, Vector avgVec, ref Vector skewVec, ref Vector kurtVec, double* end) { +#if DEBUG_ASSERT + Debug.Assert(arr + offset < end); +#endif Vector vec = VectorHelper.GetVector(arr + offset); - vec -= avgVec; + vec -= avgVec; Vector tmp = vec * vec * vec; - skewVec += tmp; - kurtVec += tmp * vec; + skewVec += tmp; + kurtVec += tmp * vec; } } } diff --git a/source/gfoidl.Stochastics/Statistics/Sample.ZTransformation.cs b/source/gfoidl.Stochastics/Statistics/Sample.ZTransformation.cs index bbd5e97..52537e1 100644 --- a/source/gfoidl.Stochastics/Statistics/Sample.ZTransformation.cs +++ b/source/gfoidl.Stochastics/Statistics/Sample.ZTransformation.cs @@ -4,6 +4,10 @@ using System.Numerics; using System.Threading.Tasks; +#if DEBUG_ASSERT +using System.Diagnostics; +#endif + namespace gfoidl.Stochastics.Statistics { partial class Sample @@ -75,61 +79,74 @@ private unsafe void ZTransformationToArrayImpl(double[] zTrans, double sigma, in fixed (double* pSource = _values) fixed (double* pTarget = zTrans) { - double* sourceArr = pSource + i; - double* targetArr = pTarget + i; + double* source = pSource + i; + double* target = pTarget + i; + n -= i; + i = 0; + double* end = source + n; if (Vector.IsHardwareAccelerated && (n - i) >= Vector.Count) { var avgVec = new Vector(avg); var sigmaInvVec = new Vector(sigmaInv); - for (; i < n - 8 * Vector.Count; i += 8 * Vector.Count) + // https://github.com/gfoidl/Stochastics/issues/46 + int m = n & ~(8 * Vector.Count - 1); + for (; i < m; i += 8 * Vector.Count) { - Core(sourceArr, targetArr, 0 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 1 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 2 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 3 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 4 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 5 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 6 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 7 * Vector.Count, avgVec, sigmaInvVec); - - sourceArr += 8 * Vector.Count; - targetArr += 8 * Vector.Count; + Core(source, target, 0 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 1 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 2 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 3 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 4 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 5 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 6 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 7 * Vector.Count, avgVec, sigmaInvVec); + + source += 8 * Vector.Count; + target += 8 * Vector.Count; } - if (i < n - 4 * Vector.Count) + m = n & ~(4 * Vector.Count - 1); + if (i < m) { - Core(sourceArr, targetArr, 0 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 1 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 2 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 3 * Vector.Count, avgVec, sigmaInvVec); - - sourceArr += 4 * Vector.Count; - targetArr += 4 * Vector.Count; - i += 4 * Vector.Count; + Core(source, target, 0 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 1 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 2 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 3 * Vector.Count, avgVec, sigmaInvVec); + + source += 4 * Vector.Count; + target += 4 * Vector.Count; + i += 4 * Vector.Count; } - if (i < n - 2 * Vector.Count) + m = n & ~(2 * Vector.Count - 1); + if (i < m) { - Core(sourceArr, targetArr, 0 * Vector.Count, avgVec, sigmaInvVec); - Core(sourceArr, targetArr, 1 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 0 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 1 * Vector.Count, avgVec, sigmaInvVec); - sourceArr += 2 * Vector.Count; - targetArr += 2 * Vector.Count; - i += 2 * Vector.Count; + source += 2 * Vector.Count; + target += 2 * Vector.Count; + i += 2 * Vector.Count; } - if (i < n - Vector.Count) + m = n & ~(1 * Vector.Count - 1); + if (i < m) { - Core(sourceArr, targetArr, 0 * Vector.Count, avgVec, sigmaInvVec); + Core(source, target, 0 * Vector.Count, avgVec, sigmaInvVec); - i += Vector.Count; + source += 1 * Vector.Count; + target += 1 * Vector.Count; } } - for (; i < n; ++i) - pTarget[i] = this.ZTransformation(pSource[i], avg, sigmaInv); + while (source < end) + { + *target = this.ZTransformation(*source, avg, sigmaInv); + source++; + target++; + } } //----------------------------------------------------------------- void Core(double* sourceArr, double* targetArr, int offset, Vector avgVec, Vector sigmaInvVec) diff --git a/source/gfoidl.Stochastics/gfoidl.Stochastics.csproj b/source/gfoidl.Stochastics/gfoidl.Stochastics.csproj index 1a45f6c..c695654 100644 --- a/source/gfoidl.Stochastics/gfoidl.Stochastics.csproj +++ b/source/gfoidl.Stochastics/gfoidl.Stochastics.csproj @@ -15,10 +15,14 @@ stochastics;statistics;outliers;error-function;Chauvenet - + bin\Release\netstandard2.0\gfoidl.Stochastics.xml + + DEBUG_ASSERT + + diff --git a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs index 96f4434..f32777c 100644 --- a/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs +++ b/tests/gfoidl.Stochastics.Tests/Statistics/SampleTests/MinMax.cs @@ -5,6 +5,7 @@ namespace gfoidl.Stochastics.Tests.Statistics.SampleTests { + [TestFixture(3)] [TestFixture(10)] [TestFixture(100)] [TestFixture(101)] @@ -36,8 +37,11 @@ public void Values___correct_MinMax() var sut = new Sample(values); - Assert.AreEqual(0, sut.Min); - Assert.AreEqual(10, sut.Max); + Assert.Multiple(() => + { + Assert.AreEqual(0, sut.Min); + Assert.AreEqual(10, sut.Max); + }); } //--------------------------------------------------------------------- [Test] @@ -54,8 +58,11 @@ public void Simd_and_ParallelizedSimd_produce_same_result() sut.GetMinMaxSimd(out double sMin, out double sMax); sut.GetMinMaxParallelizedSimd(out double pMin, out double pMax); - Assert.AreEqual(sMin, pMin, 1e-10); - Assert.AreEqual(sMax, pMax, 1e-10); + Assert.Multiple(() => + { + Assert.AreEqual(sMin, pMin, 1e-10); + Assert.AreEqual(sMax, pMax, 1e-10); + }); } //--------------------------------------------------------------------- [Test] @@ -73,8 +80,11 @@ public void Simd_and_Linq_produce_same_result() double min = values.Min(); double max = values.Max(); - Assert.AreEqual(min, actualMin, 1e-10); - Assert.AreEqual(max, actualMax, 1e-10); + Assert.Multiple(() => + { + Assert.AreEqual(min, actualMin, 1e-10); + Assert.AreEqual(max, actualMax, 1e-10); + }); } } }