diff --git a/.gitignore b/.gitignore index 28750fdb..43ac56a3 100644 --- a/.gitignore +++ b/.gitignore @@ -399,3 +399,9 @@ FodyWeavers.xsd # Cline Memory Bank - exclude from git memory-bank/ + +# macOS resource forks and metadata +._* +.DS_Store +.AppleDouble +.LSOverride diff --git a/.sonarlint/._QuanTAlib.json b/.sonarlint/._QuanTAlib.json deleted file mode 100644 index 62309917..00000000 Binary files a/.sonarlint/._QuanTAlib.json and /dev/null differ diff --git a/Directory.Build.props b/Directory.Build.props index cfdf1c70..76501f5f 100644 --- a/Directory.Build.props +++ b/Directory.Build.props @@ -44,7 +44,7 @@ - S1944,S2053,S2222,S2259,S2583,S2589,S3329,S3655,S3900,S3949,S3966,S4158,S4347,S5773,S6781 + S1944,S2053,S2222,S2259,S2583,S2589,S3329,S3655,S3776,S3900,S3949,S3966,S4158,S4347,S5773,S6781 diff --git a/lib/._QuanTAlib.Tests.csproj b/lib/._QuanTAlib.Tests.csproj deleted file mode 100644 index 62309917..00000000 Binary files a/lib/._QuanTAlib.Tests.csproj and /dev/null differ diff --git a/lib/averages/ema/Ema.Notebook.dib b/lib/averages/ema/Ema.Notebook.dib index 77bed0ba..bcf534e1 100644 --- a/lib/averages/ema/Ema.Notebook.dib +++ b/lib/averages/ema/Ema.Notebook.dib @@ -216,3 +216,79 @@ for (int i = 0; i < periods.Length; i++) } } Console.WriteLine($"\nAll Vectorized Stream/Batch values match: {allMatch}"); + +#!markdown + +## 5. Handling Invalid Values (NaN/Infinity) + +Both `Ema` and `EmaVector` use **last-value substitution** for invalid inputs. When a non-finite value (NaN, PositiveInfinity, NegativeInfinity) is encountered, it is replaced with the last valid value. This provides output continuity instead of propagating invalid values through the calculation. + +#!csharp + +Console.WriteLine("\n--- Handling Invalid Values ---"); + +// Single EMA +var emaNaN = new Ema(10); + +// Feed valid values first +emaNaN.Update(new TValue(DateTime.Now, 100.0)); +emaNaN.Update(new TValue(DateTime.Now.AddMinutes(1), 110.0)); +Console.WriteLine($"After valid values: {emaNaN.Value.Value:F2}"); + +// Feed NaN - should use last valid value (110) +var resultAfterNaN = emaNaN.Update(new TValue(DateTime.Now.AddMinutes(2), double.NaN)); +Console.WriteLine($"After NaN input: {resultAfterNaN.Value:F2} (IsFinite: {double.IsFinite(resultAfterNaN.Value)})"); + +// Feed Infinity - should use last valid value (110) +var resultAfterInf = emaNaN.Update(new TValue(DateTime.Now.AddMinutes(3), double.PositiveInfinity)); +Console.WriteLine($"After Infinity input: {resultAfterInf.Value:F2} (IsFinite: {double.IsFinite(resultAfterInf.Value)})"); + +// Continue with valid value +var resultAfterValid = emaNaN.Update(new TValue(DateTime.Now.AddMinutes(4), 120.0)); +Console.WriteLine($"After valid value (120): {resultAfterValid.Value:F2}"); + +#!csharp + +Console.WriteLine("\n--- Batch Processing with Invalid Values ---"); + +// Create series with NaN values interspersed +var seriesWithNaN = new TSeries(); +seriesWithNaN.Add(DateTime.Now.Ticks, 100.0); +seriesWithNaN.Add(DateTime.Now.Ticks + 1, 110.0); +seriesWithNaN.Add(DateTime.Now.Ticks + 2, double.NaN); +seriesWithNaN.Add(DateTime.Now.Ticks + 3, 120.0); +seriesWithNaN.Add(DateTime.Now.Ticks + 4, double.PositiveInfinity); +seriesWithNaN.Add(DateTime.Now.Ticks + 5, 130.0); + +var emaBatchNaN = new Ema(3); +var resultsWithNaN = emaBatchNaN.Update(seriesWithNaN); + +Console.WriteLine("Input → Output:"); +for (int i = 0; i < seriesWithNaN.Count; i++) +{ + var input = seriesWithNaN[i].Value; + var output = resultsWithNaN[i].Value; + var inputStr = double.IsFinite(input) ? input.ToString("F2") : input.ToString(); + Console.WriteLine($" {inputStr,-10} → {output:F2} (IsFinite: {double.IsFinite(output)})"); +} + +#!csharp + +Console.WriteLine("\n--- Vectorized EMA with Invalid Values ---"); + +int[] periodsNaN = { 5, 10 }; +var emaVectorNaN = new EmaVector(periodsNaN); + +// Feed values including invalid ones +var inputsNaN = new double[] { 100, 110, double.NaN, 120, double.PositiveInfinity, 130 }; +var time = DateTime.Now; + +foreach (var val in inputsNaN) +{ + var results = emaVectorNaN.Update(new TValue(time, val)); + var inputStr = double.IsFinite(val) ? val.ToString("F2") : val.ToString(); + Console.WriteLine($"Input: {inputStr,-10} → EMA(5): {results[0].Value:F2}, EMA(10): {results[1].Value:F2}"); + time = time.AddMinutes(1); +} + +Console.WriteLine("\nAll outputs are finite - invalid inputs were substituted with last valid values."); diff --git a/lib/averages/ema/Ema.Tests.cs b/lib/averages/ema/Ema.Tests.cs index d5dce2fd..073b7318 100644 --- a/lib/averages/ema/Ema.Tests.cs +++ b/lib/averages/ema/Ema.Tests.cs @@ -225,4 +225,106 @@ public class EmaTests Assert.Equal(100.0, result, 1e-10); } + + [Fact] + public void Ema_NaN_Input_UsesLastValidValue() + { + var ema = new Ema(10); + + // Feed some valid values + ema.Update(new TValue(DateTime.Now, 100)); + ema.Update(new TValue(DateTime.Now, 110)); + double valueBeforeNaN = ema.Value; + + // Feed NaN - should use last valid value (110) + var resultAfterNaN = ema.Update(new TValue(DateTime.Now, double.NaN)); + + // Result should be finite (not NaN) + Assert.True(double.IsFinite(resultAfterNaN.Value)); + // EMA should continue to evolve (may differ slightly due to substitution) + Assert.NotEqual(0, resultAfterNaN.Value); + } + + [Fact] + public void Ema_Infinity_Input_UsesLastValidValue() + { + var ema = new Ema(10); + + // Feed some valid values + ema.Update(new TValue(DateTime.Now, 100)); + ema.Update(new TValue(DateTime.Now, 110)); + + // Feed positive infinity - should use last valid value + var resultAfterPosInf = ema.Update(new TValue(DateTime.Now, double.PositiveInfinity)); + Assert.True(double.IsFinite(resultAfterPosInf.Value)); + + // Feed negative infinity - should use last valid value + var resultAfterNegInf = ema.Update(new TValue(DateTime.Now, double.NegativeInfinity)); + Assert.True(double.IsFinite(resultAfterNegInf.Value)); + } + + [Fact] + public void Ema_MultipleNaN_ContinuesWithLastValid() + { + var ema = new Ema(10); + + // Feed valid values + ema.Update(new TValue(DateTime.Now, 100)); + ema.Update(new TValue(DateTime.Now, 110)); + ema.Update(new TValue(DateTime.Now, 120)); + + // Feed multiple NaN values + var r1 = ema.Update(new TValue(DateTime.Now, double.NaN)); + var r2 = ema.Update(new TValue(DateTime.Now, double.NaN)); + var r3 = ema.Update(new TValue(DateTime.Now, double.NaN)); + + // All results should be finite + Assert.True(double.IsFinite(r1.Value)); + Assert.True(double.IsFinite(r2.Value)); + Assert.True(double.IsFinite(r3.Value)); + + // EMA should converge toward last valid value (120) with repeated substitution + // Values should be getting closer to 120 + Assert.True(r3.Value > r1.Value || Math.Abs(r3.Value - 120) < Math.Abs(r1.Value - 120)); + } + + [Fact] + public void Ema_BatchCalc_HandlesNaN() + { + var ema = new Ema(10); + + // Create series with NaN values interspersed + var series = new TSeries(); + series.Add(DateTime.Now.Ticks, 100); + series.Add(DateTime.Now.Ticks + 1, 110); + series.Add(DateTime.Now.Ticks + 2, double.NaN); + series.Add(DateTime.Now.Ticks + 3, 120); + series.Add(DateTime.Now.Ticks + 4, double.PositiveInfinity); + series.Add(DateTime.Now.Ticks + 5, 130); + + var results = ema.Update(series); + + // All results should be finite + foreach (var result in results) + { + Assert.True(double.IsFinite(result.Value), $"Expected finite value but got {result.Value}"); + } + } + + [Fact] + public void Ema_Reset_ClearsLastValidValue() + { + var ema = new Ema(10); + + // Feed values including NaN + ema.Update(new TValue(DateTime.Now, 100)); + ema.Update(new TValue(DateTime.Now, double.NaN)); + + // Reset + ema.Reset(); + + // After reset, first valid value should establish new baseline + var result = ema.Update(new TValue(DateTime.Now, 50)); + Assert.Equal(50.0, result.Value, 1e-10); + } } diff --git a/lib/averages/ema/Ema.cs b/lib/averages/ema/Ema.cs index b6858f49..25ba6906 100644 --- a/lib/averages/ema/Ema.cs +++ b/lib/averages/ema/Ema.cs @@ -6,9 +6,9 @@ namespace QuanTAlib; public struct EmaState { - public double Ema; - public double E; - public bool IsHot; + public double Ema { get; set; } + public double E { get; set; } + public bool IsHot { get; set; } public static EmaState New() => new() { Ema = 0, E = 1.0, IsHot = false }; } @@ -26,6 +26,7 @@ public class Ema private readonly double _alpha; private EmaState _state = EmaState.New(); private EmaState _p_state = EmaState.New(); + private double _lastValidValue; /// /// Creates EMA with specified period. @@ -36,7 +37,7 @@ public class Ema { if (period <= 0) throw new ArgumentException("Period must be greater than 0", nameof(period)); - + _alpha = 2.0 / (period + 1); } @@ -48,7 +49,7 @@ public class Ema { if (alpha <= 0 || alpha > 1) throw new ArgumentException("Alpha must be between 0 and 1", nameof(alpha)); - + _alpha = alpha; } @@ -62,22 +63,42 @@ public class Ema /// public bool IsHot => _state.IsHot; + /// + /// Gets a valid input value, using last-value substitution for non-finite inputs. + /// + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private double GetValidValue(double input) + { + if (double.IsFinite(input)) + { + _lastValidValue = input; + return input; + } + return _lastValidValue; + } + /// /// Core EMA calculation kernel. + /// Assumes input has already been validated via GetValidValue(). /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static double Compute(double input, double alpha, ref EmaState state) { state.Ema += alpha * (input - state.Ema); + double result; if (!state.IsHot) { state.E *= (1.0 - alpha); state.IsHot = state.E <= 1e-10; - return state.Ema / (1.0 - state.E); + result = state.Ema / (1.0 - state.E); } - - return state.Ema; + else + { + result = state.Ema; + } + + return result; } /// @@ -98,7 +119,9 @@ public class Ema _state = _p_state; } - double val = Compute(input.Value, _alpha, ref _state); + // Last-value substitution: replace non-finite inputs with last valid value + double val = GetValidValue(input.Value); + val = Compute(val, _alpha, ref _state); Value = new TValue(input.Time, val); return Value; } @@ -123,10 +146,12 @@ public class Ema // Local state for batch processing EmaState state = _state; - + for (int i = 0; i < len; i++) { - double val = Compute(sourceValues[i], _alpha, ref state); + // Last-value substitution: replace non-finite inputs with last valid value + double val = GetValidValue(sourceValues[i]); + val = Compute(val, _alpha, ref state); tSpan[i] = sourceTimes[i]; vSpan[i] = val; } @@ -158,6 +183,7 @@ public class Ema { _state = EmaState.New(); _p_state = _state; + _lastValidValue = 0; Value = default; } } diff --git a/lib/averages/ema/Ema.md b/lib/averages/ema/Ema.md index 7b914b05..e3b11ddd 100644 --- a/lib/averages/ema/Ema.md +++ b/lib/averages/ema/Ema.md @@ -103,6 +103,35 @@ TSeries source = ...; TSeries[] seriesResults = emaVector.Calculate(source); ``` +### Handling Invalid Values (NaN/Infinity) + +Both `Ema` and `EmaVector` use **last-value substitution** for handling invalid inputs: + +```csharp +var ema = new Ema(10); + +// Valid values establish baseline +ema.Update(new TValue(time, 100)); +ema.Update(new TValue(time, 110)); + +// NaN or Infinity inputs are replaced with last valid value (110) +var result = ema.Update(new TValue(time, double.NaN)); +Console.WriteLine(double.IsFinite(result.Value)); // true + +// Works identically for batch operations +var series = new TSeries(); +series.Add(time, 100); +series.Add(time + 1, double.NaN); // Will use 100 +series.Add(time + 2, 120); +var results = ema.Update(series); // All values are finite +``` + +**Behavior:** +- When `NaN`, `PositiveInfinity`, or `NegativeInfinity` is encountered, the last valid value is substituted +- This provides output continuity instead of propagating invalid values +- Both scalar (`Ema`) and SIMD (`EmaVector`) implementations use identical logic +- `Reset()` clears the last valid value, so the next valid input establishes a new baseline + ### Performance Characteristics * **O(1) Complexity:** The calculation time is constant regardless of the period length. diff --git a/lib/averages/ema/EmaVector.Tests.cs b/lib/averages/ema/EmaVector.Tests.cs index 5bb2b024..d0aec46f 100644 --- a/lib/averages/ema/EmaVector.Tests.cs +++ b/lib/averages/ema/EmaVector.Tests.cs @@ -135,4 +135,146 @@ public class EmaVectorTests Assert.Equal(200.0, res[0].Value, 1e-9); } + + [Fact] + public void Update_NaN_Input_UsesLastValidValue() + { + int[] periods = { 10, 20 }; + var emaVector = new EmaVector(periods); + + // Feed some valid values + emaVector.Update(new TValue(DateTime.Now, 100.0)); + emaVector.Update(new TValue(DateTime.Now, 110.0)); + + // Feed NaN - should use last valid value (110) + var resultAfterNaN = emaVector.Update(new TValue(DateTime.Now, double.NaN)); + + // All results should be finite (not NaN) + foreach (var result in resultAfterNaN) + { + Assert.True(double.IsFinite(result.Value), $"Expected finite value but got {result.Value}"); + } + } + + [Fact] + public void Update_Infinity_Input_UsesLastValidValue() + { + int[] periods = { 10, 20 }; + var emaVector = new EmaVector(periods); + + // Feed some valid values + emaVector.Update(new TValue(DateTime.Now, 100.0)); + emaVector.Update(new TValue(DateTime.Now, 110.0)); + + // Feed positive infinity - should use last valid value + var resultAfterPosInf = emaVector.Update(new TValue(DateTime.Now, double.PositiveInfinity)); + foreach (var result in resultAfterPosInf) + { + Assert.True(double.IsFinite(result.Value)); + } + + // Feed negative infinity - should use last valid value + var resultAfterNegInf = emaVector.Update(new TValue(DateTime.Now, double.NegativeInfinity)); + foreach (var result in resultAfterNegInf) + { + Assert.True(double.IsFinite(result.Value)); + } + } + + [Fact] + public void Update_MultipleNaN_ContinuesWithLastValid() + { + int[] periods = { 5, 10 }; + var emaVector = new EmaVector(periods); + + // Feed valid values + emaVector.Update(new TValue(DateTime.Now, 100.0)); + emaVector.Update(new TValue(DateTime.Now, 110.0)); + emaVector.Update(new TValue(DateTime.Now, 120.0)); + + // Feed multiple NaN values + var r1 = emaVector.Update(new TValue(DateTime.Now, double.NaN)); + var r2 = emaVector.Update(new TValue(DateTime.Now, double.NaN)); + var r3 = emaVector.Update(new TValue(DateTime.Now, double.NaN)); + + // All results should be finite + foreach (var result in r1) Assert.True(double.IsFinite(result.Value)); + foreach (var result in r2) Assert.True(double.IsFinite(result.Value)); + foreach (var result in r3) Assert.True(double.IsFinite(result.Value)); + } + + [Fact] + public void Calculate_Series_HandlesNaN() + { + int[] periods = { 5, 10 }; + var emaVector = new EmaVector(periods); + + // Create series with NaN values interspersed + var t = new System.Collections.Generic.List(); + var v = new System.Collections.Generic.List(); + var now = DateTime.Now; + + t.Add(now.Ticks); v.Add(100.0); + t.Add(now.AddMinutes(1).Ticks); v.Add(110.0); + t.Add(now.AddMinutes(2).Ticks); v.Add(double.NaN); + t.Add(now.AddMinutes(3).Ticks); v.Add(120.0); + t.Add(now.AddMinutes(4).Ticks); v.Add(double.PositiveInfinity); + t.Add(now.AddMinutes(5).Ticks); v.Add(130.0); + + var series = new TSeries(t, v); + var results = emaVector.Calculate(series); + + // All results should be finite for all periods + foreach (var periodResults in results) + { + foreach (var val in periodResults.Values) + { + Assert.True(double.IsFinite(val), $"Expected finite value but got {val}"); + } + } + } + + [Fact] + public void Reset_ClearsLastValidValue() + { + int[] periods = { 10 }; + var emaVector = new EmaVector(periods); + + // Feed values including NaN + emaVector.Update(new TValue(DateTime.Now, 100.0)); + emaVector.Update(new TValue(DateTime.Now, double.NaN)); + + // Reset + emaVector.Reset(); + + // After reset, first valid value should establish new baseline + var result = emaVector.Update(new TValue(DateTime.Now, 50.0)); + Assert.Equal(50.0, result[0].Value, 1e-9); + } + + [Fact] + public void NaN_Handling_MatchesSingleEma() + { + int[] periods = { 5, 10, 20 }; + var emaVector = new EmaVector(periods); + var emaSingles = periods.Select(p => new Ema(p)).ToArray(); + + // Data with NaN values + var values = new double[] { 10, 20, double.NaN, 40, double.PositiveInfinity, 60, 70 }; + var time = DateTime.Now; + + foreach (var val in values) + { + var tVal = new TValue(time, val); + var multiRes = emaVector.Update(tVal); + + for (int i = 0; i < periods.Length; i++) + { + var singleRes = emaSingles[i].Update(tVal); + Assert.Equal(singleRes.Value, multiRes[i].Value, 1e-9); + } + + time = time.AddMinutes(1); + } + } } diff --git a/lib/averages/ema/EmaVector.cs b/lib/averages/ema/EmaVector.cs index 54ac784e..c646fb64 100644 --- a/lib/averages/ema/EmaVector.cs +++ b/lib/averages/ema/EmaVector.cs @@ -9,6 +9,7 @@ namespace QuanTAlib; /// /// Multi-Alpha Exponential Moving Average (EMA) - SIMD optimized. /// Calculates multiple EMAs with different periods/alphas for the same input series in parallel. +/// Uses last-value substitution for invalid inputs (NaN/Infinity). /// public class EmaVector { @@ -16,7 +17,8 @@ public class EmaVector private readonly double[] _emas; private readonly double[] _Es; private readonly int _count; - + private double _lastValidValue; + /// /// Current EMA values for all periods. /// @@ -33,7 +35,7 @@ public class EmaVector _emas = new double[_count]; _Es = new double[_count]; Values = new TValue[_count]; - + for (int i = 0; i < _count; i++) { if (periods[i] <= 0) throw new ArgumentException("Period must be greater than 0", nameof(periods)); @@ -53,7 +55,7 @@ public class EmaVector _emas = new double[_count]; _Es = new double[_count]; Values = new TValue[_count]; - + for (int i = 0; i < _count; i++) { if (alphas[i] <= 0 || alphas[i] > 1) throw new ArgumentException("Alpha must be between 0 and 1", nameof(alphas)); @@ -68,6 +70,20 @@ public class EmaVector _Es[index] = 1.0; } + /// + /// Gets a valid input value, using last-value substitution for non-finite inputs. + /// + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private double GetValidValue(double input) + { + if (double.IsFinite(input)) + { + _lastValidValue = input; + return input; + } + return _lastValidValue; + } + /// /// Resets all EMA states. /// @@ -77,23 +93,27 @@ public class EmaVector { ResetAt(i); } + _lastValidValue = 0; Array.Clear(Values); } /// /// Updates EMAs with the given value. + /// Uses last-value substitution: invalid inputs (NaN/Infinity) are replaced with + /// the last known good value, providing continuity in the output series. /// /// Input value /// Array of compensated EMA values [MethodImpl(MethodImplOptions.AggressiveInlining)] public TValue[] Update(TValue input) { - double val = input.Value; - + // Last-value substitution: replace non-finite inputs with last valid value + double val = GetValidValue(input.Value); + // SIMD Loop int vecCount = Vector.Count; int i = 0; - + if (Vector.IsHardwareAccelerated && _count >= vecCount) { var vecInput = new Vector(val); @@ -110,7 +130,7 @@ public class EmaVector // Update EMA // ema += alpha * (input - ema) vecEma += vecAlpha * (vecInput - vecEma); - + // Update E (warmup factor) // E *= (1 - alpha) vecE *= (vecOne - vecAlpha); @@ -123,7 +143,7 @@ public class EmaVector // Vector.LessThanOrEqual returns Vector with all-1s for true, all-0s for false // We reinterpret as Vector for use with ConditionalSelect var isHotMask = Vector.LessThanOrEqual(vecE, vecEpsilon); - + // Select result: if hot (E <= epsilon), use raw EMA; otherwise use compensated // ConditionalSelect: mask=true -> first arg, mask=false -> second arg var vecResult = Vector.ConditionalSelect( @@ -149,7 +169,7 @@ public class EmaVector { double alpha = _alphas[i]; _emas[i] += alpha * (val - _emas[i]); - + double result = _emas[i]; if (_Es[i] > 1e-10) { @@ -159,7 +179,7 @@ public class EmaVector result = _emas[i] / (1.0 - _Es[i]); } } - + Values[i] = new TValue(input.Time, result); } @@ -175,7 +195,7 @@ public class EmaVector { int len = source.Count; var resultSeries = new TSeries[_count]; - + // Pre-allocate lists var tLists = new List[_count]; var vLists = new List[_count]; @@ -190,7 +210,7 @@ public class EmaVector var sourceValues = source.Values; var sourceTimes = source.Times; - + int vecCount = Vector.Count; var vecOne = Vector.One; var vecEpsilon = new Vector(1e-10); @@ -199,6 +219,10 @@ public class EmaVector { double val = sourceValues[t]; long time = sourceTimes[t]; + + // Last-value substitution: replace non-finite inputs with last valid value + val = GetValidValue(val); + var vecInput = new Vector(val); int i = 0; @@ -212,12 +236,12 @@ public class EmaVector vecEma += vecAlpha * (vecInput - vecEma); vecE *= (vecOne - vecAlpha); - + var vecCompensated = vecEma / (vecOne - vecE); - + // Check warmup condition: E <= 1e-10 means "hot" (use raw EMA) var isHotMask = Vector.LessThanOrEqual(vecE, vecEpsilon); - + // Select result: if hot, use raw EMA; otherwise use compensated var vecResult = Vector.ConditionalSelect( Vector.AsVectorDouble(isHotMask), @@ -241,7 +265,7 @@ public class EmaVector { double alpha = _alphas[i]; _emas[i] += alpha * (val - _emas[i]); - + double result = _emas[i]; if (_Es[i] > 1e-10) { @@ -251,7 +275,7 @@ public class EmaVector result = _emas[i] / (1.0 - _Es[i]); } } - + CollectionsMarshal.AsSpan(tLists[i])[t] = time; CollectionsMarshal.AsSpan(vLists[i])[t] = result; } diff --git a/lib/core/simd/SimdExtensions.cs b/lib/core/simd/SimdExtensions.cs index c294a1e9..66139fbd 100644 --- a/lib/core/simd/SimdExtensions.cs +++ b/lib/core/simd/SimdExtensions.cs @@ -9,15 +9,64 @@ namespace QuanTAlib; /// public static class SimdExtensions { + /// + /// Checks if span contains any non-finite values (NaN or Infinity). + /// Returns true if any non-finite value is found. + /// + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static bool ContainsNonFinite(this ReadOnlySpan span) + { + if (span.IsEmpty) return false; + + if (Vector.IsHardwareAccelerated && span.Length >= Vector.Count) + { + int vectorSize = Vector.Count; + int i = 0; + + for (; i <= span.Length - vectorSize; i += vectorSize) + { + var vector = new Vector(span.Slice(i, vectorSize)); + // Check for NaN: NaN != NaN is true + // Check for Infinity: IsInfinity + for (int j = 0; j < vectorSize; j++) + { + if (!double.IsFinite(vector[j])) + return true; + } + } + + // Check remaining elements + for (; i < span.Length; i++) + { + if (!double.IsFinite(span[i])) + return true; + } + + return false; + } + + // Scalar fallback + for (int i = 0; i < span.Length; i++) + { + if (!double.IsFinite(span[i])) + return true; + } + return false; + } + /// /// Calculates sum using SIMD vectorization when available. /// 4-8x faster than scalar loop on AVX2/AVX-512 hardware. + /// Returns NaN if any input value is non-finite. /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static double SumSIMD(this ReadOnlySpan span) { if (span.IsEmpty) return 0.0; + // Guard against non-finite inputs + if (span.ContainsNonFinite()) return double.NaN; + if (Vector.IsHardwareAccelerated && span.Length >= Vector.Count) { Vector sum = Vector.Zero; @@ -53,6 +102,7 @@ public static class SimdExtensions /// /// Calculates minimum value using SIMD vectorization when available. /// 4-6x faster than scalar loop on AVX2/AVX-512 hardware. + /// Returns NaN if any input value is non-finite. /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static double MinSIMD(this ReadOnlySpan span) @@ -60,6 +110,9 @@ public static class SimdExtensions if (span.IsEmpty) return double.NaN; if (span.Length == 1) return span[0]; + // Guard against non-finite inputs + if (span.ContainsNonFinite()) return double.NaN; + if (Vector.IsHardwareAccelerated && span.Length >= Vector.Count) { int vectorSize = Vector.Count; @@ -104,6 +157,7 @@ public static class SimdExtensions /// /// Calculates maximum value using SIMD vectorization when available. /// 4-6x faster than scalar loop on AVX2/AVX-512 hardware. + /// Returns NaN if any input value is non-finite. /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static double MaxSIMD(this ReadOnlySpan span) @@ -111,6 +165,9 @@ public static class SimdExtensions if (span.IsEmpty) return double.NaN; if (span.Length == 1) return span[0]; + // Guard against non-finite inputs + if (span.ContainsNonFinite()) return double.NaN; + if (Vector.IsHardwareAccelerated && span.Length >= Vector.Count) { int vectorSize = Vector.Count; @@ -155,25 +212,34 @@ public static class SimdExtensions /// /// Calculates average using SIMD vectorization when available. /// 4-8x faster than scalar loop on AVX2/AVX-512 hardware. + /// Returns NaN if any input value is non-finite. /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static double AverageSIMD(this ReadOnlySpan span) { if (span.IsEmpty) return double.NaN; + // SumSIMD already guards against non-finite, which will propagate NaN return span.SumSIMD() / span.Length; } /// /// Calculates variance using SIMD vectorization (Welford's online algorithm adapted). /// More numerically stable than naive two-pass algorithm. + /// Returns NaN if any input value is non-finite or if mean is non-finite. /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static double VarianceSIMD(this ReadOnlySpan span, double? mean = null) { if (span.Length < 2) return double.NaN; + // Guard against non-finite inputs + if (span.ContainsNonFinite()) return double.NaN; + double m = mean ?? span.AverageSIMD(); + // Guard against non-finite mean (could be passed in or computed from non-finite values) + if (!double.IsFinite(m)) return double.NaN; + if (Vector.IsHardwareAccelerated && span.Length >= Vector.Count) { var meanVec = new Vector(m); @@ -216,16 +282,19 @@ public static class SimdExtensions /// /// Calculates standard deviation using SIMD vectorization. + /// Returns NaN if any input value is non-finite or if mean is non-finite. /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static double StdDevSIMD(this ReadOnlySpan span, double? mean = null) { + // VarianceSIMD already guards against non-finite, which will propagate NaN through Sqrt return Math.Sqrt(span.VarianceSIMD(mean)); } /// /// Finds both min and max in a single pass using SIMD vectorization. /// More efficient than calling MinSIMD and MaxSIMD separately. + /// Returns (NaN, NaN) if any input value is non-finite. /// [MethodImpl(MethodImplOptions.AggressiveInlining)] public static (double Min, double Max) MinMaxSIMD(this ReadOnlySpan span) @@ -233,6 +302,9 @@ public static class SimdExtensions if (span.IsEmpty) return (double.NaN, double.NaN); if (span.Length == 1) return (span[0], span[0]); + // Guard against non-finite inputs + if (span.ContainsNonFinite()) return (double.NaN, double.NaN); + if (Vector.IsHardwareAccelerated && span.Length >= Vector.Count) { int vectorSize = Vector.Count;