在x64架构中更快的sin()函数

4

主要问题

有没有人有一个适用于x64架构的快速sin()实现?它不需要是纯Pascal。

说明

我有一个VCL应用程序,有些情况下当它为x64编译时运行得比较慢。

它执行了大量的浮点3D计算,并且我已经跟踪到这个问题是由于当输入值变大时,System.Sin()System.Cos()在x64上变得非常缓慢。

我通过创建一个简单的测试应用程序来计时,该应用程序测量计算sin(x)所需的时间,使用不同的x值进行测试,发现其性能差异巨大:

              call:     x64:     x86:
              Sin(1)   16 ms    20 ms
             Sin(10)   30 ms    20 ms
            Sin(100)   32 ms    20 ms
           Sin(1000)   34 ms    21 ms
          Sin(10000)   30 ms    21 ms
         Sin(100000)   30 ms    16 ms
        Sin(1000000)   35 ms    20 ms
       Sin(10000000)  581 ms    20 ms
      Sin(100000000) 1026 ms    21 ms
     Sin(1000000000) 1187 ms    22 ms
    Sin(10000000000) 1320 ms    21 ms
   Sin(100000000000) 1456 ms    20 ms
  Sin(1000000000000) 1581 ms    17 ms
 Sin(10000000000000) 1717 ms    22 ms
Sin(100000000000000) 1846 ms    23 ms
           Sin(1E15) 1981 ms    21 ms
           Sin(1E16) 2100 ms    21 ms
           Sin(1E17) 2240 ms    22 ms
           Sin(1E18) 2372 ms    18 ms
                etc    etc      etc

您在这里看到的是sin(1E5)运行速度大约是sin(1E8)的300倍。

如果您感兴趣,我是这样创建上述表格的:

{$APPTYPE CONSOLE}
program SinTest;

uses Diagnostics, Math, SysUtils;

var
  i : Integer;
  x : double;
  sw: TStopwatch;

begin
  x := 1;

  while X < 1E18 do
  begin
    sw    := TStopwatch.StartNew;
    for i := 1 to 500000 do
      System.Sin(x);

    // WriteLn(System.sin(x), #9,System.Sin(fmod(x,2*pi)));

    sw.Stop;

    WriteLn('    ', ('Sin(' + round(x).ToString + ')'):20, ' ', sw.ElapsedMilliseconds,' ms');

    x := x * 10;
  end;

  WriteLn('Press any key to continue');
  readln;
end.

注意事项:

  • 在StackOverflow上有一些关于更快的正弦函数的问题,但它们中没有一个有用于转换为Delphi的源代码,就像这个:C++中最快的正弦、余弦和平方根实现(不需要太高的准确度)

  • x64的其余部分比32位的运行得更快。

  • 我发现了一些垃圾解决方法,就是这样做:Sin(FMod(x,2*pi))。它可以提供正确的结果,并且对于较大的数字运行速度很快。当然,对于较小的数字,它会慢一些。


2
你可能并不关心精度,否则你就不会使用如此大的值来调用三角函数了。你肯定知道舍入误差意味着对于这样的输入值,三角函数是没有意义的吧?或者对你来说准确性并不重要? - David Heffernan
1
那么,看看你能否猜出这个程序的输出:{$APPTYPE CONSOLE} var s1, s2: Single; begin s1 := 10000000.5; s2 := 10000000.0; Writeln(s1=s2); end. 这里有一个提示。输出不是 FALSE - David Heffernan
1
看起来 MSVC 可以更快地完成它,我很想知道是如何做到的,因为我敢打赌对于合理的输入值,它也可以更快地完成。但对于您的大型输入值,即使调用这些三角函数也是浪费时间的,正如我之前的评论所示。 - David Heffernan
2
不,你正在使用单精度。这在问题中已经有了。或者你的问题不是你想问的那个吗? - David Heffernan
1
@LURD:也许我应该将这样的东西封装成一个新的sin()函数,然后使用它。 - Wouter van Nifterick
显示剩余14条评论
2个回答

3

虽然在用户模式代码中这可能会被强烈反对(在内核模式代码中完全禁止),但如果您希望在x64代码中保留传统的x87行为,您可以编写以下函数:

function SinX87(x:double):double;
var
  d : double;
asm
  movsd qword ptr [rbp+8], xmm0
  fld qword ptr [rbp+8]
  fsin
  fstp qword ptr [rbp+8]
  movsd xmm0, qword ptr [rbp+8]
end;

这会增加一些开销,因为您需要将该值从SSE寄存器弹出到堆栈中,将其加载到x87单元中进行计算,将该值弹回到堆栈中,然后将其重新加载到XMM0以获得函数结果。不过,sin的计算非常繁重,因此这只是一个相对较小的开销。我只会在需要保留x87的sin实现的任何特殊性时才这样做。

存在其他库可以比Delphi的纯Pascal例程更有效地计算x64代码中的sin。我的首选是导出一组良好的C ++例程到DLL中。 此外,正如David所说,使用具有荒谬大参数的三角函数也不是一个明智的选择。


很酷,速度非常稳定,无论输入是什么。对于小于π的值,它稍微慢一点;其他情况下总是更快。结果与Delphi的System.Sin()略有不同,但对我需要处理的数字来说并不重要。结果看起来很好。这正是我所需要的。现在我只需要添加一些丑陋的{$ifdef}代码,就可以恢复在x64下的性能了。谢谢! - Wouter van Nifterick
@WoutervanNifterick 另外,我不确定异常会如何处理... 我肯定会先进行测试。在 x64 模式下,x87 控制字是否会默认设置为任何合理的值也不确定 - 我很快就完成了这个,但要注意一些警告。 - J...
测试过了,确实处理方式有些不同。例如SinX87(NaN)不会像System.Sin()一样引发任何异常。因此确实存在差异,但这是一个很好的帮助。我将进行一些额外的测试,但到目前为止,它看起来完全符合我的需求。 - Wouter van Nifterick

2
如果你对我的最终解决方案感兴趣:
我进行了一些实验,通过这种方式(如LU RD和e所述)。Jerry Coffin建议的方法:
function sin(x:double):double;
begin
  if x<1E6 then
    Result := system.sin(x)
  else
    Result := system.sin(fmod(x,2*pi));
end;

也许这与我的特定CPU上测试代码的可预测性有关,但如果我不使用if,只是始终使用fmod(),则计算较小的值实际上会更快。奇怪的是,一些除法需要进行,这应该比比较两个值要慢才对。
所以现在我最终使用的是:
function sin(const x: double): double; { inline; }
begin
  {$IFDEF CPUX64}
  Result := System.sin(Math.FMod(x,2*pi));
  {$ELSE}
  Result := System.sin(x);
  {$ENDIF}
end;

顺便提一下,添加inline后,它的运行速度甚至快了1.5倍。在我的机器上,它的速度与J...函数完全相同。但即使没有使用Inline,它已经比System.Sin()快了数百倍,所以我选择使用它。


1
即使你使用fmod(x, 2*pi),如@DavidHeffernan所指出的那样,你也会遇到一个事实,即作为双精度变量的x不能容纳超过约17位小数位的信息,因此你失去了传递到sin函数中的精度。例如:如果你将x从100000000000000000.0步进到100000000000000000.1,表示0.1弧度步长,这两个数字是相同的,因为当0.1被添加时,它会丢失,因为双精度变量不足以容纳整个值。你必须找到另一种编码x的方法。 - Mike Dunlavey

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接