午後わてんのブログ

ベランダ菜園とWindows用アプリ作成(WPFとC#)

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

f:id:gogowaten:20200212114500p:plain

{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

f:id:gogowaten:20200218150749p:plain

要素の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

 

 

f:id:gogowaten:20200218155040p:plain

素数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

 

 計測結果をグラフで

f:id:gogowaten:20200218155006p:plain

一覧

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型にした。

 

 

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

f:id:gogowaten:20200218171957p:plain

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

gogowaten.hatenablog.com

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

f:id:gogowaten:20200218172122p:plain

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

 

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

f:id:gogowaten:20200218172444p:plain

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

 

Test02は偏差の引き算をint型

f:id:gogowaten:20200218173005p:plain

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

 

System.Numerics.Vector.Dotメソッド

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

f:id:gogowaten:20200212115205p:plain

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

f:id:gogowaten:20200212120659p:plain
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メソッドを使ってみると

f:id:gogowaten:20200212121248p:plain

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

 

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

Vector<float>を使ったTest03

f:id:gogowaten:20200218173551p:plain

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

f:id:gogowaten:20200219121629p:plain

今の環境で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>で計算してオーバーフロー

f:id:gogowaten:20200218174655p:plain

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

 

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

f:id:gogowaten:20200218175412p:plain

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

f:id:gogowaten:20200218200251p:plain

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

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

 

Test12

f:id:gogowaten:20200218201026p:plain

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

 

Vectorでその2

f:id:gogowaten:20200218201824p:plain

Test03をその2にしただけ

 

Vector4

f:id:gogowaten:20200218205144p:plain

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

 

Vector.Widen

f:id:gogowaten:20200218203441p:plain

拡大変換

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

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

とかできる

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

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

f:id:gogowaten:20200218203154p:plain

変換が多いから遅くなりそうだったけど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