StdDev.cs

using System.Runtime.Intrinsics;
using System.Runtime.Intrinsics.X86;
using BenchmarkDotNet.Attributes;
using BenchmarkDotNet.Jobs;
using BenchmarkDotNet.Running;

#pragma warning disable CS8618

var summary = BenchmarkRunner.Run<MainTest>();

[SimpleJob(RuntimeMoniker.Net70, baseline: true)]
[WarmupCount(2)]
[IterationCount(2)]
[MemoryDiagnoser]
public class MainTest
{
    const int ITERATIONS = 100;

    [Params(64, 4096, 32768)]
    public int N { get; set; }

    double[] src;

    [IterationSetup]
    public void Setup()
    {
        src = Enumerable.Range(0, N)
                        .Select(p => Random.Shared.NextDouble() * (p + 1))
                        .ToArray();
    }

    [Benchmark(Baseline = true)]
    public double StdDev()
    {
        var r = 0.0;
        for (var iteration = 0; iteration < ITERATIONS; iteration++)
        {
            r = 0;

            var mean     = src.Average();
            var variance = src.Select(p => Math.Pow(p - mean, 2)).Average();
            r += Math.Sqrt(variance);
        }

        return r;
    }

    [Benchmark]
    public unsafe double StdDevVector()
    {
        var r = 0.0;
        for (var iteration = 0; iteration < ITERATIONS; iteration++)
        {
            fixed (double* srcPtr = src)
            {
                var pointer = srcPtr;

                var vsum = Vector256.Load(pointer);
                pointer += 4;
                for (var i = 0; i < src.Length / 4; i++, pointer += 4)
                    vsum = Avx.Add(vsum, Vector256.Load(pointer));

                var vectorMean = Vector256.Create((vsum.GetElement(0) +
                                                   vsum.GetElement(1) +
                                                   vsum.GetElement(2) +
                                                   vsum.GetElement(3)) / src.Length);

                pointer = srcPtr;
                var acc = 0.0;
                for (var i = 0; i < src.Length / 4; i++, pointer += 4)
                {
                    var diff = Avx.Subtract(Vector256.Load(pointer), vectorMean);
                    var pow  = Avx.Multiply(diff, diff);
                    acc += pow.GetElement(0) +
                           pow.GetElement(1) +
                           pow.GetElement(2) +
                           pow.GetElement(3);
                }

                r += Math.Sqrt(acc / src.Length);
            }
        }

        return r;
    }
}