C#.NETで分散を求める2回め、System.Numerics.Vector.Dotも使ってみたけど…

 今回の記事は前回の

gogowaten.hatenablog.com

↑この記事を書き直したもので

直した理由は速度比較をしていたんだけど、Math.Powの使い所が間違っていて、速度比較になっていないことが記事を書いた6日後に

gogowaten.hatenablog.com

わかったから、なので前回のコードからMathPowを一掃して最初から計測し直しただけなので、内容的にはほとんど同じ

分散の求め方は二通りあるようで

mathtrain.jp

こちらを参考にした

分散その1 = 偏差の2乗の平均(偏差は平均との差のこと)

分散その2 = 2乗の平均 - 平均の2乗

分散その1

{20, 21, 7, 12}

っていう4個の要素のデータ、これの分散は33.5

エクセル2007だとVARPっていう関数で求めることができる

分散 = 偏差の2乗の平均値

偏差 = 各要素と平均値との差

分散 = 各要素と平均値との差の2乗の平均値

データの平均値

{20, 21, 7, 12}の平均値

(20+21+7+12) / 4

= 15

偏差

{20, 21, 7, 12}の各偏差は(各要素 - 平均値)なので

{20-15, 21-15, 7-15, 12-15}

= {5, 6, -8, -3}

偏差の2乗

{5*5, 6*6, -8*-8, -3*-3}

= {25, 36, 64, 9}

分散は偏差の2乗の平均値

(25+36+64+9) / 4

= 33.5

分散その2

要素の2乗の平均は、(20*20+21*21+7*7+12*12)/4=258.5

平均の2乗は、平均が(20+21+7+12)/4=15、これを2乗して15*15=225

要素の2乗の平均 - 平均の2乗は、258.5-225=33.5

要素数1千万のbyte型配列の分散を100回求める処理時間計測

配列の値はアプリ起動時にランダムなのを入れるので、分散の結果は起動ごとに変わる

前回の記事との違いは

• Test01とTest11で使っていたMathPowをやめて、普通の掛け算にした
• それによって速くなったので配列の要素数を10倍にした、100万→1千万
• 一斉テストボタンを追加した

アプリの制作と計測環境

• CPU AMD Ryzen 5 2400G(4コア8スレッド)
• MEM DDR4-2666
• Window 10 Home 64bit
• Visual Studio 2019 Community .NET Core 3.1 WPF C#

.NET Frameworkだと参照の追加とかあってめんどくさいけど、.NET Coreならラク

アプリのダウンロード先とコードはギットハブ

ファイル名：20200210_配列の分散、VectorDot_v1.1.zip

github.com

 計測結果をグラフで

一覧

Test	分散	精度	Vector型	速度
Test01	その1	小数点		1.00	Vectorを使わない普通の計算をしているこれ基準にして比較
Test02	その1	整数		1.27	最速、整数計算だから精度が落ちるけど要素数が多いので誤差
Test03	その1	小数点	float	0.42	偏差とドット積をVectorを使って計算、Vector自体は速いはずだけど準備に時間かかって遅い
Test04	その1	小数点	double	0.26	↑とほぼ同じ
Test05	その1	整数	int	0.51	↑とほぼ同じ
Test06	その1	破綻	byte	0.95	偏差とトッド積の計算でオーバーフローなので求めることができない
Test07	その1	破綻	short	0.49	↑とほぼ同じ
Test08	その1	小数点	Vector4	0.90	Vectorのなかでは速いけど普通の計算のほうが速い
Test09	その1	整数	int	0.58	偏差は普通に計算してドット積だけVector.Dot
Test10	その1	小数点	double	0.28	↑とほぼ同じ
Test11	その2	小数点		1.09	ここからその2で計算、Test01の改変だけど少し速くなった
Test12	その2	小数点		1.07	Test02をその2にした、遅くなったけど精度は上
Test13	その2	小数点	float	0.45	Test03をその2にした、速度は0.42から微増
Test14	その2	小数点	double	0.28	Test04をその2にした、速度は0.26から微増
Test15	その2	小数点	int	0.52	Test05をその2にした、速度は0.51から微増
Test16	その2	破綻	byte	0.97	Test06をその2にした、速度は0.95から微増
Test17	その2	破綻	short	0.50	Test07をその2にした、速度は0.49から微増
Test18	その2	小数点	Vector4	0.96	Test08をその2にした、速度は0.90から微増
Test19	その2	小数点	uint	0.85	VectorのWidenメソッドでVector<byte>をushort→uintにして計算Vectorのなかでは速いけど
Test20	その2	小数点	uint	0.87	↑の変数の使い方を少し変えただけ、意味なかった

残念ながらVectorを使わないほうが速いという結果になった、けどAVXが強いIntelのCPUやZen2コアになってAVXが強化されたAMD RyzenならVectorのほうが速いかも？

使うとしたら最速のTest02、精度が気になるならTest11かTest12かな、Vectorの出番はない。

分散の求め方その1とその2の比較では、ほんの少しだけどその2のほうが速い

要素数1千万、これは写真撮影に使っているスマカメの画素数が800万だから、これに近い切りのいい数字を選んだ。画像処理で使う予定だから配列の型もbyte型にした。

平均値は先に計算しておいた

ParallelForEachとPartitionerを使ったマルチスレッドは前回の方法

gogowaten.hatenablog.com

アプリ起動時にこれで計算して

MyAveargeって名前つけたフィールドに置いておいて使うことにした

速度の基準にするTest01は普通の掛け算

前回はべき乗計算でMathPowを使っていたけど、普通の掛け算にした、これだけで30倍以上速くなった、これのせいで記事を書き直す羽目になった

Test02は偏差の引き算をint型

double型で計算していたTest01をint型で計算するようにしたTest02、精度を少し犠牲にして少し速度アップ。前回もここではMathPowを使わないほうが速いことに気づいていたんだけど、このときはint型で計算したために速くなったと思っていた

System.Numerics.Vector.Dotメソッド

System.NumericsはCPUのSIMDを使って計算するクラスがいっぱい

ドット積はググってみたけどよくわからん

2つの同じ要素数の配列(ベクトル？)の各値を掛け算して、それを合計することみたい。掛け算して足し算

配列1{1, 3, 5, 7}

配列2{16, 8, 2, 4}

この2つの配列のドット積は78

各要素同士の掛け算

{1*16, 3*8, 5*2, 7*4} = {16,24,10,28}

合計

{16+24+10+28} = 78

ドット積は78

Dotメソッドを使ってみると

期待通りの値78が返ってきた

分散は偏差の2乗の平均、このうち偏差の2乗は掛け算だし、平均を求めるときには合計する。つまり掛け算して足し算、これはドット積(Dotメソッド)と同じ！

Vector<float>を使ったTest03

偏差(平均値との差)はVector.Subtract、ドット積はVector.Dotを使った計算

今の環境でVector<float>.Countは8、1回で8個同時に計算できるから、Vectorの計算自体は速いはずだけど、121行目から124行目でVector用の配列を作るのに時間がかかって全体では遅くなるみたいで、処理速度は基準の0.42倍、2倍以上時間がかかる

このTest03、04、05までは型が違うだけで同じ傾向

あと、今回もVectorでの計算は少し問題があって、VectorはVectorCountの値の要素数ごとにしか計算できないから、全要素数がVectorCountで割り切れなかった場合、この余りは別に計算しないと正しい値にならないんだけど、めんどくさいので今回は計算していない例えば全要素数が10000でVectorCountが6だと、10000%6=4で、この4個の余りは計算していない

Test06はVector<byte>で計算してオーバーフロー

もとの配列の型と同じbyte型だから、新たに配列を作成する必要がない分速いけど、byte型にマイナスの値がないから偏差でオーバーフロー、ドット積でもオーバーフローで間違った値になる。配列から直接Vectorを作ったら速いんじゃないかなあと思ったんだけど、速度は基準の0.95倍、Vector<byte>は32個同時に計算だからもっと速くなってもいいと思うんだけどねえ

Vector.Dotは同じ型で返すからすぐオーバーフローする

byte型VectorでのDotはbyte型で返す、今の環境ではbyte型Vectorは一度に32個を計算するから、かんたんにbyte型の最大値255を超えてしまう

偏差が全部たったの3だったとしてもドット積は288になってしまう、3*3=9、これが32個だから、9*32=288

これがあるからもとの配列の型がbyte型でも、longやdoubleの配列に変換してVector<long>やVector<double>にして計算する必要がある、その変換が時間かかってトータルではVectorを使わないほうが速いってなるんだけど、使い方あるのかしらねえ

Test11、ここからは分散の求め方その2

求め方その1だったTest01をその2にしたのがこのTest11、前回記事ではここでもMathPowを使っていたので今回は掛け算に変更したことで30倍以上速くなった。299行目だったのを300行目に変更

使うならこのTest11と次のTest12、速くて精度も高いし、書く文字数も少なくて済む

Test12

Test02をその2にしただけ、ing型での計算だけどbyte型はもともと整数だから、2乗和の集計の変数の型はlongのほうが良かったかもって今思った、これはTest11でもそうだなあ

Vectorでその2

Test03をその2にしただけ

Vector4

一度に計算できる要素数は4個、型はfloatと決まっているけどVectorの中では最速

Vector.Widen

拡大変換

1つのVector<byte>を2つのVector<ushort>に変換

1つのVector<ushort>なら2つのVector<uint>に変換

とかできる

できないのは符号付きから符号なしの型やその逆、整数型から小数点型やその逆とかはできない、残念

それでもbyteからuintまで変換できればオーバーフローしなくなるので、これを使ったのがTest19と20

変換が多いから遅くなりそうだったけどVectorを使った中では速いほう

<Window x:Class="_20200210_配列の分散_VectorDot.MainWindow"
        xmlns="http://schemas.microsoft.com/winfx/2006/xaml/presentation"
        xmlns:x="http://schemas.microsoft.com/winfx/2006/xaml"
        xmlns:d="http://schemas.microsoft.com/expression/blend/2008"
        xmlns:mc="http://schemas.openxmlformats.org/markup-compatibility/2006"
        xmlns:local="clr-namespace:_20200210_配列の分散_VectorDot"
        mc:Ignorable="d"
        Title="MainWindow" Height="600" Width="614">
  <Grid>
    <StackPanel>
      <StackPanel.Resources>
        <Style TargetType="StackPanel">
          <Setter Property="Margin" Value="2"/>
        </Style>
        <Style TargetType="Button">
          <Setter Property="Width" Value="60"/>
        </Style>
        <Style TargetType="TextBlock">
          <Setter Property="Margin" Value="2,0"/>
        </Style>
      </StackPanel.Resources>
      <TextBlock x:Name="MyTextBlock" Text="text" HorizontalAlignment="Center"/>
      <TextBlock x:Name="MyTextBlockVectorCount" Text="vectorCount" HorizontalAlignment="Center"/>
      <TextBlock x:Name="MyTextBlockCpuThreadCount" Text="threadCount" HorizontalAlignment="Center"/>
      <StackPanel Orientation="Horizontal" HorizontalAlignment="Center">
        <Button x:Name="ButtonAll" Content="一斉テスト" Margin="20,0" Width="120"/>
        <TextBlock x:Name="TbAll" Text="time"/>
        <Button x:Name="ButtonReset" Content="reset" Margin="20,0"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button1" Content="test1"/>
        <TextBlock x:Name="Tb1" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button2" Content="test2"/>
        <TextBlock x:Name="Tb2" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button3" Content="test3"/>
        <TextBlock x:Name="Tb3" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button4" Content="test4"/>
        <TextBlock x:Name="Tb4" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button5" Content="test5"/>
        <TextBlock x:Name="Tb5" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button6" Content="test6"/>
        <TextBlock x:Name="Tb6" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button7" Content="test7"/>
        <TextBlock x:Name="Tb7" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button8" Content="test8"/>
        <TextBlock x:Name="Tb8" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button9" Content="test9"/>
        <TextBlock x:Name="Tb9" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button10" Content="test10"/>
        <TextBlock x:Name="Tb10" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button11" Content="test11"/>
        <TextBlock x:Name="Tb11" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button12" Content="test12"/>
        <TextBlock x:Name="Tb12" Text="time"/>
      </StackPanel>
      <!--<Border Height="1" Background="Orange" UseLayoutRounding="True"/>-->
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button13" Content="test13"/>
        <TextBlock x:Name="Tb13" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button14" Content="test14"/>
        <TextBlock x:Name="Tb14" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button15" Content="test15"/>
        <TextBlock x:Name="Tb15" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button16" Content="test16"/>
        <TextBlock x:Name="Tb16" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button17" Content="test17"/>
        <TextBlock x:Name="Tb17" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button18" Content="test18"/>
        <TextBlock x:Name="Tb18" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button19" Content="test19"/>
        <TextBlock x:Name="Tb19" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button20" Content="test20"/>
        <TextBlock x:Name="Tb20" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button21" Content="test21"/>
        <TextBlock x:Name="Tb21" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button22" Content="test22"/>
        <TextBlock x:Name="Tb22" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button23" Content="test23"/>
        <TextBlock x:Name="Tb23" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button24" Content="test24"/>
        <TextBlock x:Name="Tb24" Text="time"/>
      </StackPanel>
      <StackPanel Orientation="Horizontal">
        <Button x:Name="Button25" Content="test25"/>
        <TextBlock x:Name="Tb25" Text="time"/>
      </StackPanel>



    </StackPanel>
  </Grid>
</Window>

using System;
using System.Threading.Tasks;
using System.Windows;
using System.Windows.Controls;
using System.Diagnostics;
using System.Numerics;
using System.Threading;
using System.Collections.Concurrent;
using System.Collections.Generic;


namespace _20200210_配列の分散_VectorDot
{
    /// <summary>
    /// Interaction logic for MainWindow.xaml
    /// </summary>
    public partial class MainWindow : Window
    {
        private byte[] MyByteAry;
        //private int[] MyIntAry;
        //private long[] MyLongAry;
        private const int LOOP_COUNT = 100;
        private const int ELEMENT_COUNT = 10_000_000;
        private double MyAverage;

        public MainWindow()
        {
            InitializeComponent();



            MyTextBlock.Text = $"byte型配列の値の分散、要素数{ELEMENT_COUNT.ToString("N0")}の分散を{LOOP_COUNT}回求める";
            MyTextBlockVectorCount.Text = $"Vector<long>.Count = {Vector<long>.Count}";
            string str = $"VectorCount : Long={Vector<long>.Count}, Double={Vector<double>.Count}, int={Vector<int>.Count}, flort={Vector<float>.Count}, short={Vector<short>.Count}, byte={Vector<byte>.Count}";
            MyTextBlockVectorCount.Text = str;
            MyTextBlockCpuThreadCount.Text = $"CPUスレッド数：{Environment.ProcessorCount.ToString()} thread";
            //MyIntAry = Enumerable.Range(1, ELEMENT_COUNT).ToArray();//連番値
            //MyIntAry = Enumerable.Repeat(1, ELEMENT_COUNT).ToArray();//全値1
            //MyLongAry = new long[MyIntAry.Length];//long型配列作成
            //MyIntAry.CopyTo(MyLongAry, 0);
            MyInitialize();

            ButtonAll.Click += async (s, e) => await MyExeAll();
            ButtonReset.Click += (s, e) => MyReset();

            Button1.Click += (s, e) => MyExe(Test01_Double_ForLoop, Tb1, MyByteAry);
            Button2.Click += (s, e) => MyExe(Test02_Integer_ForLoop, Tb2, MyByteAry);
            Button3.Click += (s, e) => MyExe(Test03_FloatVectorSubtractDot, Tb3, MyByteAry);
            Button4.Click += (s, e) => MyExe(Test04_DoubleVectorSubtractDot, Tb4, MyByteAry);
            Button5.Click += (s, e) => MyExe(Test05_IntegerVectorSubtractDot, Tb5, MyByteAry);
            Button6.Click += (s, e) => MyExe(Test06_ByteVectorSubtractDot_Overflow, Tb6, MyByteAry);
            Button7.Click += (s, e) => MyExe(Test07_ShortVectorSubtractDot_Overflow, Tb7, MyByteAry);
            Button8.Click += (s, e) => MyExe(Test08_FloatVector4, Tb8, MyByteAry);
            Button9.Click += (s, e) => MyExe(Test09_IntegerVectorDot, Tb9, MyByteAry);
            Button10.Click += (s, e) => MyExe(Test10_DoubleVectorDot, Tb10, MyByteAry);

            Button11.Click += (s, e) => MyExe(Test11_Double_ForLoop, Tb11, MyByteAry);
            Button12.Click += (s, e) => MyExe(Test12_Integer_ForLoop, Tb12, MyByteAry);
            Button13.Click += (s, e) => MyExe(Test13_FloatVectorDot, Tb13, MyByteAry);
            Button14.Click += (s, e) => MyExe(Test14_DoubleVectorDot, Tb14, MyByteAry);
            Button15.Click += (s, e) => MyExe(Test15_IntegerVectorDot, Tb15, MyByteAry);
            Button16.Click += (s, e) => MyExe(Test16_ByteVectorDot_Overflow, Tb16, MyByteAry);
            Button17.Click += (s, e) => MyExe(Test17_ShortVectorDot_Overflow, Tb17, MyByteAry);
            Button18.Click += (s, e) => MyExe(Test18_FloatVector4, Tb18, MyByteAry);

            Button19.Click += (s, e) => MyExe(Test19_Byte_ushort_uintVectorDot, Tb19, MyByteAry);
            Button20.Click += (s, e) => MyExe(Test20_Byte_ushort_uintVectorDot, Tb20, MyByteAry);
        }
        private void MyInitialize()
        {
            MyByteAry = new byte[ELEMENT_COUNT];
            var r = new Random();
            r.NextBytes(MyByteAry);
            //要素の平均値
            //MyByteAry = new byte[] { 0, 21, 7, 255 };
            MyAverage = GetAverage(MyByteAry);
        }

        private double Test01_Double_ForLoop(byte[] ary)
        {
            //平均との差(偏差)の2乗を合計
            double total = 0;
            double temp;
            for (int i = 0; i < ary.Length; i++)
            {
                //total += Math.Pow(ary[i] - MyAverage, 2.0);//遅すぎ
                temp = ary[i] - MyAverage;
                total += temp * temp;
            }
            //合計 / 要素数 = 分散
            return total / ary.Length;
        }
        private double Test02_Integer_ForLoop(byte[] ary)
        {
            //平均との差の2乗を合計
            long total = 0;
            int average = (int)MyAverage;
            int ii;//ループの外に出したほうが誤差程度に速い
            for (int i = 0; i < ary.Length; i++)
            {
                //total += (int)Math.Pow(ary[i] - average, 2);
                //total += (ary[i] - average) * (ary[i] - average);//こっちのほうが↑より10倍以上速い
                ii = ary[i] - average;
                total += ii * ii;
            }
            //合計 / 要素数 = 分散
            return total / (double)ary.Length;
        }

        //Vector<float>で計算
        private double Test03_FloatVectorSubtractDot(byte[] ary)
        {
            var vAverage = new Vector<float>((float)MyAverage);
            int simdLength = Vector<float>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<float> v;
            var ss = new float[simdLength];
            double total = 0;
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                //平均との差
                v = System.Numerics.Vector.Subtract(vAverage, new Vector<float>(ss));
                //差の2乗を合計
                total += System.Numerics.Vector.Dot(v, v);
            }
            return total / ary.Length;
        }

        //Vector<double>で計算
        private double Test04_DoubleVectorSubtractDot(byte[] ary)
        {
            var vAverage = new Vector<double>(MyAverage);
            int simdLength = Vector<double>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<double> v;
            var ss = new double[simdLength];
            double total = 0;
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                //平均との差
                v = System.Numerics.Vector.Subtract(vAverage, new Vector<double>(ss));
                //差の2乗を合計
                total += System.Numerics.Vector.Dot(v, v);
            }
            return total / ary.Length;
        }

        //Vector<int>で計算
        private double Test05_IntegerVectorSubtractDot(byte[] ary)
        {
            var vAverage = new Vector<int>((int)MyAverage);
            int simdLength = Vector<int>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<int> v;
            var ss = new int[simdLength];
            double total = 0;
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                //平均との差
                v = System.Numerics.Vector.Subtract(vAverage, new Vector<int>(ss));
                //差の2乗を合計
                total += System.Numerics.Vector.Dot(v, v);
            }
            return total / ary.Length;
        }

        //Vector<byte>で計算、ドット積でオーバーフロー
        private double Test06_ByteVectorSubtractDot_Overflow(byte[] ary)
        {
            var vAverage = new Vector<byte>((byte)MyAverage);
            int simdLength = Vector<byte>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<byte> v;
            double total = 0;

            for (int i = 0; i < lastIndex; i += simdLength)
            {
                //平均との差は、byte型にはマイナスがないから間違った値になる
                v = System.Numerics.Vector.Subtract(vAverage, new Vector<byte>(ary, i));
                //差の2乗を合計
                total += System.Numerics.Vector.Dot(v, v);//オーバーフロー
            }
            return total / ary.Length;
        }

        //Vector<short>で計算はドット積でオーバーフロー
        private double Test07_ShortVectorSubtractDot_Overflow(byte[] ary)
        {
            var vAverage = new Vector<short>((short)MyAverage);
            int simdLength = Vector<short>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<short> v;
            long total = 0;
            var ss = new short[simdLength];
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                //平均との差
                v = System.Numerics.Vector.Subtract(vAverage, new Vector<short>(ss));
                //差の2乗を合計
                total += System.Numerics.Vector.Dot(v, v);//オーバーフロー
            }
            return (double)total / ary.Length;
        }

        //Vector4
        private double Test08_FloatVector4(byte[] ary)
        {
            int lastIndex = ary.Length - (ary.Length % 4);
            var vAverage = new Vector4((float)MyAverage);
            Vector4 v;//= new Vector4();
            double total = 0;
            for (int i = 0; i < lastIndex; i += 4)
            {
                //平均との差
                v = Vector4.Subtract(new Vector4(ary[i], ary[i + 1], ary[i + 2], ary[i + 3]), vAverage);
                //差の2乗を合計
                total += Vector4.Dot(v, v);
            }
            return total / ary.Length;
        }

        //平均との差は普通のループで、int型
        private double Test09_IntegerVectorDot(byte[] ary)
        {
            int average = (int)MyAverage;
            int simdLength = Vector<int>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<int> v;
            long total = 0;
            int[] ii = new int[simdLength];
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                //平均との差を配列に入れる
                for (int j = 0; j < simdLength; j++)
                {
                    ii[j] = average - ary[i + j];
                }
                //差の配列からVector作成して2乗和を合計していく
                v = new Vector<int>(ii);
                total += System.Numerics.Vector.Dot(v, v);
            }
            return total / (double)ary.Length;
        }

        //平均との差は普通のループで、double型
        private double Test10_DoubleVectorDot(byte[] ary)
        {
            int simdLength = Vector<double>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<double> v;
            double total = 0;
            double[] ii = new double[simdLength];
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                //平均との差の配列
                for (int j = 0; j < simdLength; j++)
                {
                    ii[j] = MyAverage - ary[i + j];
                }
                //差の配列からVector作成してドット積
                v = new Vector<double>(ii);
                total += System.Numerics.Vector.Dot(v, v);
            }
            return total / ary.Length;
        }





        #region 分散 = 2乗和の平均 - 平均の2乗
        //分散の意味と求め方、分散公式の使い方
        //https://sci-pursuit.com/math/statistics/variance.html



        private double Test11_Double_ForLoop(byte[] ary)
        {
            //2乗和
            double total = 0;
            for (int i = 0; i < ary.Length; i++)
            {
                //total += Math.Pow(ary[i], 2.0);//遅すぎ
                total += ary[i] * ary[i];
            }
            //2乗和の平均
            total /= ary.Length;
            //2乗和の平均 - 平均の2乗
            return total - (MyAverage * MyAverage);// Math.Pow(MyAverage, 2.0);
        }

        private double Test12_Integer_ForLoop(byte[] ary)
        {

            double total = 0;
            int ii;
            for (int i = 0; i < ary.Length; i++)
            {
                ii = ary[i];
                total += ii * ii;
            }
            total /= ary.Length;

            return total - (MyAverage * MyAverage);
        }

        //float
        private double Test13_FloatVectorDot(byte[] ary)
        {
            int simdLength = Vector<float>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<float> v;
            var ss = new float[simdLength];
            double total = 0;
            for (int i = 0; i < lastIndex; i += simdLength)
            //配列を作成してVector作成してドット積
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                v = new Vector<float>(ss);
                total += System.Numerics.Vector.Dot(v, v);
            }
            //2乗和の平均 - 平均の2乗
            return (total / ary.Length) - (MyAverage * MyAverage);
        }


        //Vector<double>で計算
        private double Test14_DoubleVectorDot(byte[] ary)
        {
            int simdLength = Vector<double>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<double> v;
            var ss = new double[simdLength];
            double total = 0;
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                v = new Vector<double>(ss);
                total += System.Numerics.Vector.Dot(v, v);
            }
            return (total / ary.Length) - (MyAverage * MyAverage);
        }

        //Vector<int>で計算
        private double Test15_IntegerVectorDot(byte[] ary)
        {
            int simdLength = Vector<int>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<int> v;
            var ss = new int[simdLength];
            double total = 0;
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                v = new Vector<int>(ss);
                total += System.Numerics.Vector.Dot(v, v);
            }
            return (total / ary.Length) - (MyAverage * MyAverage);
        }

        //Vector<byte>で計算、ドット積でオーバーフロー
        private double Test16_ByteVectorDot_Overflow(byte[] ary)
        {
            int simdLength = Vector<byte>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<byte> v;
            double total = 0;

            for (int i = 0; i < lastIndex; i += simdLength)
            {
                v = new Vector<byte>(ary, i);
                total += System.Numerics.Vector.Dot(v, v);//オーバーフロー
            }
            return (total / ary.Length) - (MyAverage * MyAverage);
        }

        //Vector<short>で計算はドット積でオーバーフロー
        private double Test17_ShortVectorDot_Overflow(byte[] ary)
        {
            int simdLength = Vector<short>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<short> v;
            long total = 0;
            var ss = new short[simdLength];
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                for (int j = 0; j < simdLength; j++)
                {
                    ss[j] = ary[i + j];
                }
                v = new Vector<short>(ss);
                total += System.Numerics.Vector.Dot(v, v);//オーバーフロー
            }
            return (total / ary.Length) - (MyAverage * MyAverage);
        }

        //Vector4
        private double Test18_FloatVector4(byte[] ary)
        {
            int lastIndex = ary.Length - (ary.Length % 4); ;
            Vector4 v;
            double total = 0;
            for (int i = 0; i < lastIndex; i += 4)
            {
                v = new Vector4(ary[i], ary[i + 1], ary[i + 2], ary[i + 3]);
                total += Vector4.Dot(v, v);
            }
            return (total / ary.Length) - (MyAverage * MyAverage);
        }



        //Vector<byte>をWidenでVector<ushort>にしてドット積計算はオーバーフローだったので
        //Vector<ushort>からさらにVector<uint>にしてドット積
        private double Test19_Byte_ushort_uintVectorDot(byte[] ary)
        {
            int simdLength = Vector<byte>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<byte> v;
            double total = 0;
            Vector<ushort> v1; Vector<ushort> v2;
            Vector<uint> vv1; Vector<uint> vv2; Vector<uint> vv3; Vector<uint> vv4;
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                v = new Vector<byte>(ary, i);
                System.Numerics.Vector.Widen(v, out v1, out v2);
                System.Numerics.Vector.Widen(v1, out vv1, out vv2);
                System.Numerics.Vector.Widen(v2, out vv3, out vv4);
                total += System.Numerics.Vector.Dot(vv1, vv1);
                total += System.Numerics.Vector.Dot(vv2, vv2);
                total += System.Numerics.Vector.Dot(vv3, vv3);
                total += System.Numerics.Vector.Dot(vv4, vv4);
            }
            return (total / ary.Length) - (MyAverage * MyAverage);
        }

        //↑と同じ、インライン化しただけ
        private double Test20_Byte_ushort_uintVectorDot(byte[] ary)
        {
            int simdLength = Vector<byte>.Count;
            int lastIndex = ary.Length - (ary.Length % simdLength);
            Vector<byte> v;
            double total = 0;
            for (int i = 0; i < lastIndex; i += simdLength)
            {
                v = new Vector<byte>(ary, i);
                System.Numerics.Vector.Widen(v, out Vector<ushort> v1, out Vector<ushort> v2);
                System.Numerics.Vector.Widen(v1, out Vector<uint> vv1, out Vector<uint> vv2);
                System.Numerics.Vector.Widen(v2, out Vector<uint> vv3, out Vector<uint> vv4);
                total += System.Numerics.Vector.Dot(vv1, vv1);
                total += System.Numerics.Vector.Dot(vv2, vv2);
                total += System.Numerics.Vector.Dot(vv3, vv3);
                total += System.Numerics.Vector.Dot(vv4, vv4);
            }
            return (total / ary.Length) - (MyAverage * MyAverage);
        }


        #endregion




        //一斉テスト
        private async Task MyExeAll()
        {
            this.IsEnabled = false;
            var sw = new Stopwatch();
            sw.Start();
            await Task.Run(() => MyExe(Test01_Double_ForLoop, Tb1, MyByteAry));
            await Task.Run(() => MyExe(Test02_Integer_ForLoop, Tb2, MyByteAry));
            await Task.Run(() => MyExe(Test03_FloatVectorSubtractDot, Tb3, MyByteAry));
            await Task.Run(() => MyExe(Test04_DoubleVectorSubtractDot, Tb4, MyByteAry));
            await Task.Run(() => MyExe(Test05_IntegerVectorSubtractDot, Tb5, MyByteAry));
            await Task.Run(() => MyExe(Test06_ByteVectorSubtractDot_Overflow, Tb6, MyByteAry));
            await Task.Run(() => MyExe(Test07_ShortVectorSubtractDot_Overflow, Tb7, MyByteAry));
            await Task.Run(() => MyExe(Test08_FloatVector4, Tb8, MyByteAry));
            await Task.Run(() => MyExe(Test09_IntegerVectorDot, Tb9, MyByteAry));
            await Task.Run(() => MyExe(Test10_DoubleVectorDot, Tb10, MyByteAry));
            await Task.Run(() => MyExe(Test11_Double_ForLoop, Tb11, MyByteAry));
            await Task.Run(() => MyExe(Test12_Integer_ForLoop, Tb12, MyByteAry));
            await Task.Run(() => MyExe(Test13_FloatVectorDot, Tb13, MyByteAry));
            await Task.Run(() => MyExe(Test14_DoubleVectorDot, Tb14, MyByteAry));
            await Task.Run(() => MyExe(Test15_IntegerVectorDot, Tb15, MyByteAry));
            await Task.Run(() => MyExe(Test16_ByteVectorDot_Overflow, Tb16, MyByteAry));
            await Task.Run(() => MyExe(Test17_ShortVectorDot_Overflow, Tb17, MyByteAry));
            await Task.Run(() => MyExe(Test18_FloatVector4, Tb18, MyByteAry));
            await Task.Run(() => MyExe(Test19_Byte_ushort_uintVectorDot, Tb19, MyByteAry));
            await Task.Run(() => MyExe(Test20_Byte_ushort_uintVectorDot, Tb20, MyByteAry));

            sw.Stop();
            this.IsEnabled = true;
            TbAll.Text = $"一斉テスト処理時間：{sw.Elapsed.TotalSeconds.ToString("000.000")}秒";
        }

        //全てのtextblockのTextを空白にする
        private void MyReset()
        {
            var tbs = GetObjects<TextBlock>(this);
            foreach (var item in tbs)
            {
                item.Text = "";
            }
        }
        //コントロールの列挙
        private static List<T> GetObjects<T>(DependencyObject obj) where T : DependencyObject
        {
            var lo = new List<T>();
            foreach (var item in LogicalTreeHelper.GetChildren(obj))
            {
                if (item is T) { lo.Add((T)item); }
                if(item is DependencyObject dObj)
                {
                    foreach (var item2 in GetObjects<T>(dObj))
                    {
                        lo.Add(item2);
                    }
                }
            }
            return lo;
        }


        //平均値
        private double GetAverage(byte[] ary)
        {
            long total = 0;
            Parallel.ForEach(Partitioner.Create(0, ary.Length),
                (range) =>
                {
                    long subtotal = 0;
                    for (int i = range.Item1; i < range.Item2; i++)
                    {
                        subtotal += ary[i];
                    }
                    Interlocked.Add(ref total, subtotal);
                });
            return total / (double)ary.Length;
        }










        private void MyExe(Func<byte[], double> func, TextBlock tb, byte[] ary)
        {
            double total = 0;
            var sw = new Stopwatch();
            sw.Start();
            for (int i = 0; i < LOOP_COUNT; i++)
            {
                total = func(ary);
            }
            sw.Stop();
            this.Dispatcher.Invoke(() =>
            {
                tb.Text = $"処理時間：{sw.Elapsed.TotalSeconds.ToString("00.000")}秒  分散 = {total.ToString("F4")}  {System.Reflection.RuntimeReflectionExtensions.GetMethodInfo(func).Name}";
            });
        }



    }

}

 次回はマルチスレッドで比較したい

関連記事

次は2日後

8日前

gogowaten.hatenablog.com