ドット積を求めるSystem.Numerics.VectorのDotメソッドを使って分散を求めてみた
2020/02/19追記ここから
処理時間の基準にしたTest01でMathPowの使い方に問題があることが6日後の
↑の記事でわかったので
書き直したのが ↓
なのでここから↓はほとんど意味ない
2020/02/19追記ここまで
結果
普通のループやSystem.Numerics.Vector、Vector4を使って分散を求めてみた
値はランダムなbyte型、要素数100万から100回
アプリのダウンロード先
ギットハブ
環境
CPU AMD Ryzen 5 2400G(4コア8スレッド)
MEM DDR4-2666
Window 10 Home 64bit
Visual Studio 2019 Community .NET Core 3.1 WPF C#
分散の求め方その1
{20, 21, 7, 12}
っていう4個の要素のデータ、これの分散は33.5
分散 = 偏差の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
System.Numerics.Vector.Dotメソッド
ドット積はググってみたけどよくわからんけど
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メソッドを使ってみると
期待通りの値が返ってきた
分散は偏差の2乗の平均、このうち偏差の2乗は掛け算だし、平均を求めるときには合計する。つまり掛け算して足し算、これはドット積(Dotメソッド)と同じ！
byte型配列の分散を求める
値はランダムなbyte型の配列
平均値は先に計算しておく
//平均値 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; }
型はdoubleにした
ParallelForEachとPartitionerを使ったマルチスレッドは前回の方法
平均値は
アプリの起動時に計算してフィールドに置いておく、32行目
普通のループで計算、処理時間は3.4秒
private double Test01_Double_ForLoop(byte[] ary) { //平均との差(偏差)の2乗を合計 double total = 0; for (int i = 0; i < ary.Length; i++) { total += Math.Pow(ary[i] - MyAverage, 2.0); } //合計 / 要素数 = 分散 return total / ary.Length; }
System.Numerics.Vector.Dotメソッドを使って計算、処理時間は0.3秒
//平均との差は普通のループで、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; }
少し問題があって、VectorはVectorCountの値の要素数ごとにしか計算できないから、全要素数がVectorCountで割り切れなかった場合、この余りは別に計算しないと正しい値にならないんだけど、めんどくさいので今回は計算していない
全要素数が10でVectorCountが4のとき、10%4=2で2個のあまりが出る、
System.Numerics.Vector.Subtractメソッドと
System.Numerics.Vector.Dotメソッドを使って計算、処理時間は0.3秒
//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; }
Subtractメソッドは要素同士の引き算、これで平均値との差を求めている
ここまでの3つだと普通のループよりもDotメソッドを使ったほうが10倍速い！ってところなんだけど
普通のループをint型で計算、処理時間0.07秒
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; }
平均値を切り捨ての整数にして、全てをint型で計算したら0.07秒！整数だから正確な分散からは少し違う値になるけど
ランダムなbyte型の値100万の場合
5474.5395 double
5474.7733 int
その差0.2338、5000分の0.2しか違わない、これならint型で十分
じゃあVectorもint型で計算したら
//平均との差は普通のループで、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; }
これでの処理時間は0.15秒、double型Vectorより速くなったけど、それでも普通のループで計算したほうが2倍速い
引き算もVectorで
//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; }
これでの処理時間は0.18、これも普通のループのほうが速い
int型での普通のループが速すぎてVectorの出番がない
分散の求め方その2
ここみると別の求め方があるみたいで
要素の2乗の平均 - 平均の2乗
これで分散になるみたい、説明見ても全然わからんけど
たしかに分散になる、すごい
これを使って普通のループで
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); } //2乗和の平均 total /= ary.Length; //2乗和の平均 - 平均の2乗 return total - 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); }
上がdouble型、下がint型で計算
結果はそれぞれ求め方その1からその2で
3.36秒から2.84秒
0.070秒から0.085秒
double型は少し速くなったけど、int型では少し遅くなったって、今見たら最後の計算でdouble型になってるからこれが原因かも
return total - (MyAverage * MyAverage);を
return total - ((int)MyAverage * (int)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); }
0.343から0.332秒 double型
0.154から0.187秒 int型
これもint型が少し速くなったけど最後の平均の2乗がdouble型で計算しているせいかも？
どちらにしてもあんまり変わらない結果
Widenメソッドも使って
//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); }
最初のVectorはbyte型配列から直接作成できるbyte型のVector、これをWidenメソッドでushort型Vectorに変換、これにもWidenメソッドを使ってuint型Vectorに変換してドット積
結果は0.107秒、速くなったけど、これでも普通のint型ループの0.070秒にはかなわなかった
System.Numerics.Vector4クラス
Vector4クラスはfloat型の値を4つで色々計算できるみたいで、SubtractやDotメソッドもある、これを使って普通の分散
//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; }
0.101秒！速い
Vector4で別の求め方の分散
//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); }
0.097秒！測定はある程度ブレるから0.1秒台も出るけど0.09秒台もでる、今回のテストでVectorを使った中では最速、float型なのでint型よりも正確な分散のはずなので、速さと正確さの両方を求めるならこれが良さそう
結果をグラフで
double型での普通ループの1番と11番は遅すぎるので除外して
分散の求め方その1がtest01から10
分散の求め方その2がtest11から20で概ねその2のほうが速い
overflowって付いてる赤いグラフ6,7,16,17番は、計算途中で値がオーバーフローして正確ではない値が返ってくるので、速くても意味ない参考値
Vectorの型による速度は速い順に
byte＞＞short＝int＞float＞＞doubleになった
ただしbyteとshortはオーバーフローするから意味ないし、intだと気分的に誤差が不安なのでfloatがいい、そしてVectorでは最速だったVector4はちょうどfloat型なので、使うとしたら18番かな。整数でいいって場合なら2番
それにしてもint型で計算する普通のループ(2番)がこんなに速いとは思わなかったねえ、48倍だもんなあ、それでもdouble型同士ならVectorのほうが10倍だし、十分な精度だと思うfloat型ならVector4で35倍だからSIMDを使うSystem.Numerics.Vectorは速いってことで、あと気になるのはVectorだと要素数がVector.Countで割り切れなかった場合のあまりの要素の処理は普通のループですることになるからそれはめんどくさいねえ
今回は細切れコードを画像じゃなくて文字で貼り付けてみたけど、画像を貼り付けるほうがラクかも、コピペするときは行頭に無駄なタブがあるから、それを削除して貼り付けるのがめんどくさかった
MainWindow.xaml
<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"> <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>
MainWindow.xaml.cs
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; 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 = 1_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(); 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[] { 20, 21, 7, 12 }; MyAverage = GetAverage(MyByteAry); } private double Test01_Double_ForLoop(byte[] ary) { //平均との差(偏差)の2乗を合計 double total = 0; for (int i = 0; i < ary.Length; i++) { total += Math.Pow(ary[i] - MyAverage, 2.0); } //合計 / 要素数 = 分散 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); } //2乗和の平均 total /= ary.Length; //2乗和の平均 - 平均の2乗 return total - 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 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(); tb.Text = $"処理時間：{sw.Elapsed.TotalSeconds.ToString("00.000")}秒 分散 = {total.ToString("F4")} {System.Reflection.RuntimeReflectionExtensions.GetMethodInfo(func).Name}"; } } }
関連記事
次回は6日後
前回は昨日
