C#中的数学优化

Question

我整天都在分析一个应用程序,并且已经优化了几个代码,我在todo列表中留下了这个.它是神经网络的激活函数,被调用超过1亿次.根据dotTrace,它占总功能时间的约60%.

你会如何优化这个？

public static float Sigmoid(double value) {
    return (float) (1.0 / (1.0 + Math.Pow(Math.E, -value)));
}

Answer 1

尝试:

public static float Sigmoid(double value) {
    return 1.0f / (1.0f + (float) Math.Exp(-value));
}

编辑:我做了一个快速的基准测试.在我的机器上,上面的代码比你的方法快了大约43%,这个数学上等效的代码是最快的一点(比原来快46%):

public static float Sigmoid(double value) {
    float k = Math.Exp(value);
    return k / (1.0f + k);
}

编辑2:我不确定C#函数有多少开销,但如果你#include <math.h>在源代码中,你应该能够使用它,它使用float-exp函数.它可能会快一点.

public static float Sigmoid(double value) {
    float k = expf((float) value);
    return k / (1.0f + k);
}

此外,如果您正在进行数百万次调用,则函数调用开销可能是个问题.尝试制作内联函数,看看是否有任何帮助.

Answer 2

如果它是激活函数,如果e ^ x的计算完全准确,那么它是否非常重要？

例如,如果您使用近似值(1 + x/256)^ 256,在我的Pentium Pentium测试中(我假设C#基本上编译为相同的处理器指令),这比e ^ x快大约7-8倍(Math.exp()),精确到2位小数,最大约为x +/- 1.5,并且在您所述范围内的正确数量级内.(显然,要提高到256,你实际上将数字平方8次 - 不要使用Math.Pow!)在Java中:

double eapprox = (1d + x / 256d);
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;

保持加倍或减半256(并添加/删除乘法)取决于您希望逼近的准确度.即使n = 4,对于-0.5和0.5之间的x值,它仍然给出大约1.5的小数位精度(并且看起来比Math.exp()快15倍).

PS我忘了提-你应该显然不是真正的256分:乘以一个常数1/256.Java的JIT编译器自动进行这种优化(至少Hotspot会这样做),我假设C#也必须这样做.

Answer 3

看看这篇文章.它有一个用Java编写的e ^ x的近似值,它应该是它的C#代码(未经测试):

public static double Exp(double val) {  
    long tmp = (long) (1512775 * val + 1072632447);  
    return BitConverter.Int64BitsToDouble(tmp << 32);  
}

在我的基准测试中,这比Math.exp()(在Java中)快5倍.该近似基于论文" 指数函数的快速,紧凑近似 ",其被精确地开发用于神经网络.它基本上与2048条目的查找表和条目之间的线性近似相同,但所有这些都与IEEE浮点技巧有关.

编辑:根据Special Sauce,它比CLR实现快3.25倍.谢谢!

Answer 4

1. 请记住,此激活函数的任何更改都会以不同行为为代价.这甚至包括切换到浮动(从而降低精度)或使用激活替代品.只有尝试使用您的用例才能显示正确的方法.
2. 除了简单的代码优化之外,我还建议考虑计算的并行化(即:利用机器的多个内核甚至Windows Azure云上的机器)并改进训练算法.

更新: 在ANN激活功能的查找表上发布

更新2:我删除了LUT上的点,因为我把它们与完全散列混淆了.谢谢Henrik Gustafsson让我回到赛道上.所以内存不是问题,尽管搜索空间仍然会被局部极值弄乱.

Answer 5

在1亿次呼叫时,我开始怀疑探查器开销是否会影响您的结果.用no-op替换计算,看看是否仍然报告占用了60%的执行时间......

或者更好的是,创建一些测试数据并使用秒表计时器来分析大约一百万个呼叫.

Answer 6

如果你能够与C++互操作,你可以考虑将所有值存储在一个数组中并使用SSE循环它们,如下所示:

void sigmoid_sse(float *a_Values, float *a_Output, size_t a_Size){
    __m128* l_Output = (__m128*)a_Output;
    __m128* l_Start  = (__m128*)a_Values;
    __m128* l_End    = (__m128*)(a_Values + a_Size);

    const __m128 l_One        = _mm_set_ps1(1.f);
    const __m128 l_Half       = _mm_set_ps1(1.f / 2.f);
    const __m128 l_OneOver6   = _mm_set_ps1(1.f / 6.f);
    const __m128 l_OneOver24  = _mm_set_ps1(1.f / 24.f);
    const __m128 l_OneOver120 = _mm_set_ps1(1.f / 120.f);
    const __m128 l_OneOver720 = _mm_set_ps1(1.f / 720.f);
    const __m128 l_MinOne     = _mm_set_ps1(-1.f);

    for(__m128 *i = l_Start; i < l_End; i++){
        // 1.0 / (1.0 + Math.Pow(Math.E, -value))
        // 1.0 / (1.0 + Math.Exp(-value))

        // value = *i so we need -value
        __m128 value = _mm_mul_ps(l_MinOne, *i);

        // exp expressed as inifite series 1 + x + (x ^ 2 / 2!) + (x ^ 3 / 3!) ...
        __m128 x = value;

        // result in l_Exp
        __m128 l_Exp = l_One; // = 1

        l_Exp = _mm_add_ps(l_Exp, x); // += x

        x = _mm_mul_ps(x, x); // = x ^ 2
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_Half, x)); // += (x ^ 2 * (1 / 2))

        x = _mm_mul_ps(value, x); // = x ^ 3
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver6, x)); // += (x ^ 3 * (1 / 6))

        x = _mm_mul_ps(value, x); // = x ^ 4
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver24, x)); // += (x ^ 4 * (1 / 24))

#ifdef MORE_ACCURATE

        x = _mm_mul_ps(value, x); // = x ^ 5
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver120, x)); // += (x ^ 5 * (1 / 120))

        x = _mm_mul_ps(value, x); // = x ^ 6
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver720, x)); // += (x ^ 6 * (1 / 720))

#endif

        // we've calculated exp of -i
        // now we only need to do the '1.0 / (1.0 + ...' part
        *l_Output++ = _mm_rcp_ps(_mm_add_ps(l_One,  l_Exp));
    }
}

但是,请记住,您将使用的数组应使用_aligned_malloc(some_size*sizeof(float),16)进行分配,因为SSE需要将内存与边界对齐.

使用SSE,我可以在大约半秒内计算所有1亿个元素的结果.但是,一次分配那么多内存将花费你几乎三分之一的千兆字节,所以我建议一次处理更多但更小的数组.您甚至可能想要考虑使用100K或更多元素的双缓冲方法.

此外,如果元素的数量开始大幅增长,您可能希望选择在GPU上处理这些内容(只需创建一维float4纹理并运行一个非常简单的片段着色器).

Answer 7

FWIW,这是我已经发布的答案的C#基准.(Empty是一个只返回0的函数,用于测量函数调用开销)

Empty Function:       79ms   0
Original:             1576ms 0.7202294
Simplified: (soprano) 681ms  0.7202294
Approximate: (Neil)   441ms  0.7198783
Bit Manip: (martinus) 836ms  0.72318
Taylor: (Rex Logan)   261ms  0.7202305
Lookup: (Henrik)      182ms  0.7204863

public static object[] Time(Func<double, float> f) {
    var testvalue = 0.9456;
    var sw = new Stopwatch();
    sw.Start();
    for (int i = 0; i < 1e7; i++)
        f(testvalue);
    return new object[] { sw.ElapsedMilliseconds, f(testvalue) };
}
public static void Main(string[] args) {
    Console.WriteLine("Empty:       {0,10}ms {1}", Time(Empty));
    Console.WriteLine("Original:    {0,10}ms {1}", Time(Original));
    Console.WriteLine("Simplified:  {0,10}ms {1}", Time(Simplified));
    Console.WriteLine("Approximate: {0,10}ms {1}", Time(ExpApproximation));
    Console.WriteLine("Bit Manip:   {0,10}ms {1}", Time(BitBashing));
    Console.WriteLine("Taylor:      {0,10}ms {1}", Time(TaylorExpansion));
    Console.WriteLine("Lookup:      {0,10}ms {1}", Time(LUT));
}

Answer 8

注:这是一个后续这一职位.

编辑:更新计算同样的事情,这个和这个,采取一些灵感此.

现在看看你让我做了什么!你让我安装了Mono!

$ gmcs -optimize test.cs && mono test.exe
Max deviation is 0.001663983
10^7 iterations using Sigmoid1() took 1646.613 ms
10^7 iterations using Sigmoid2() took 237.352 ms

C再也不值得付出努力了,世界正在前进:)

所以,速度提高了10 6倍.有窗口框的人可以使用MS-stuff调查内存使用情况和性能:)

使用LUT进行激活功能并不少见,特别是在硬件中实现时.如果您愿意包含这些类型的表,那么该概念有许多经过充分验证的变体.但是,正如已经指出的那样,混叠可能会成为一个问题,但也有办法解决这个问题.进一步阅读:

• NEURObjects由乔治·瓦伦蒂尼(有也是在这一个纸)
• 具有数字LUT激活功能的神经网络
• 通过降低精度激活函数来提升神经网络特征提取
• 一种新的具有整数权重和量化非线性激活函数的神经网络学习算法
• 量化对高阶函数神经网络的影响

有些人对此有所了解:

• 当你到达桌子外面时错误会增加(但在极端情况下会收敛到0); 对于x约+ -7.0.这是由于选择的缩放因子.较大的SCALE值在中间范围内给出较高的误差,但在边缘处较小.
• 这通常是一个非常愚蠢的测试,我不知道C#,它只是我的C代码的简单转换:)
• Rinat Abdullin非常正确,别名和精度损失可能会导致问题,但由于我没有看到变量,我只能建议你试试这个.事实上,我同意他所说的一切,除了查找表的问题.

请原谅复制粘贴编码......

using System;
using System.Diagnostics;

class LUTTest {
    private const float SCALE = 320.0f;
    private const int RESOLUTION = 2047;
    private const float MIN = -RESOLUTION / SCALE;
    private const float MAX = RESOLUTION / SCALE;

    private static readonly float[] lut = InitLUT();

    private static float[] InitLUT() {
      var lut = new float[RESOLUTION + 1];

      for (int i = 0; i < RESOLUTION + 1; i++) {
        lut[i] = (float)(1.0 / (1.0 + Math.Exp(-i / SCALE)));
      }
      return lut;
    }

    public static float Sigmoid1(double value) {
        return (float) (1.0 / (1.0 + Math.Exp(-value)));
    }

    public static float Sigmoid2(float value) {
      if (value <= MIN) return 0.0f;
      if (value >= MAX) return 1.0f;
      if (value >= 0) return lut[(int)(value * SCALE + 0.5f)];
      return 1.0f - lut[(int)(-value * SCALE + 0.5f)];
    }

    public static float error(float v0, float v1) {
      return Math.Abs(v1 - v0);
    }

    public static float TestError() {
        float emax = 0.0f;
        for (float x = -10.0f; x < 10.0f; x+= 0.00001f) {
          float v0 = Sigmoid1(x);
          float v1 = Sigmoid2(x);
          float e = error(v0, v1);
          if (e > emax) emax = e;
        }
        return emax;
    }

    public static double TestPerformancePlain() {
        Stopwatch sw = new Stopwatch();
        sw.Start();
        for (int i = 0; i < 10; i++) {
            for (float x = -5.0f; x < 5.0f; x+= 0.00001f) {
                Sigmoid1(x);
            }
        }
        sw.Stop();
        return sw.Elapsed.TotalMilliseconds;
    }    

    public static double TestPerformanceLUT() {
        Stopwatch sw = new Stopwatch();
        sw.Start();
        for (int i = 0; i < 10; i++) {
            for (float x = -5.0f; x < 5.0f; x+= 0.00001f) {
                Sigmoid2(x);
            }
        }
        sw.Stop();
        return sw.Elapsed.TotalMilliseconds;
    }    

    static void Main() {
        Console.WriteLine("Max deviation is {0}", TestError());
        Console.WriteLine("10^7 iterations using Sigmoid1() took {0} ms", TestPerformancePlain());
        Console.WriteLine("10^7 iterations using Sigmoid2() took {0} ms", TestPerformanceLUT());
    }
}

Answer 9

在.NET数学算法中,F#比C#具有更好的性能.因此,在F#中重写神经网络可能会提高整体性能.

如果我们在F#中重新实现LUT基准测试片段(我一直在使用稍微调整过的版本),那么结果代码:

• 在588.8ms而不是3899,2ms执行sigmoid1基准测试
• 在156.6ms而不是 411.4ms执行sigmoid2(LUT)基准测试

更多细节可以在博客文章中找到.这是F#片段JIC:

#light

let Scale = 320.0f;
let Resolution = 2047;

let Min = -single(Resolution)/Scale;
let Max = single(Resolution)/Scale;

let range step a b =
  let count = int((b-a)/step);
  seq { for i in 0 .. count -> single(i)*step + a };

let lut = [| 
  for x in 0 .. Resolution ->
    single(1.0/(1.0 +  exp(-double(x)/double(Scale))))
  |]

let sigmoid1 value = 1.0f/(1.0f + exp(-value));

let sigmoid2 v = 
  if (v <= Min) then 0.0f;
  elif (v>= Max) then 1.0f;
  else
    let f = v * Scale;
    if (v>0.0f) then lut.[int (f + 0.5f)]
    else 1.0f - lut.[int(0.5f - f)];

let getError f = 
  let test = range 0.00001f -10.0f 10.0f;
  let errors = seq { 
    for v in test -> 
      abs(sigmoid1(single(v)) - f(single(v)))
  }
  Seq.max errors;

open System.Diagnostics;

let test f = 
  let sw = Stopwatch.StartNew(); 
  let mutable m = 0.0f;
  let result = 
    for t in 1 .. 10 do
      for x in 1 .. 1000000 do
        m <- f(single(x)/100000.0f-5.0f);
  sw.Elapsed.TotalMilliseconds;

printf "Max deviation is %f\n" (getError sigmoid2)
printf "10^7 iterations using sigmoid1: %f ms\n" (test sigmoid1)
printf "10^7 iterations using sigmoid2: %f ms\n" (test sigmoid2)

let c = System.Console.ReadKey(true);

和输出(针对F#1.9.6.2 CTP的Release编译,没有调试器):

Max deviation is 0.001664
10^7 iterations using sigmoid1: 588.843700 ms
10^7 iterations using sigmoid2: 156.626700 ms

更新:更新基准测试以使用10 ^ 7次迭代使结果与C相媲美

UPDATE2:以下是来自同一台机器的C实现的性能结果:

Max deviation is 0.001664
10^7 iterations using sigmoid1: 628 ms
10^7 iterations using sigmoid2: 157 ms