Delphi Tokyo 64 位将非正规数刷新为零？

Posted 2023-02-22

【中文标题】Delphi Tokyo 64 位将非正规数刷新为零？

【英文标题】：Delphi Tokyo 64-bit flushes denormal numbers to zero?

【发布时间】：2018-07-20 15:11:58

【问题描述】：

在短暂查看 system.math 的源代码时，我发现 从下面的程序可以看出，64 位版本的 Delphi Tokyo 10.2.3 将非规范 IEEE-Double 刷新为零；

$apptype console
uses
  system.sysutils, system.math;
var
  x: double;
const
  twopm1030 : UInt64 = $0000100000000000; 2^(-1030)
begin
  x := PDouble(@twopm1030)^;
  writeln(x);
  x := ldexp(1,-515);
  writeln(x*x);
  x := ldexp(1,-1030);
  writeln(x);
end.

对于 32 位，输出符合预期

8.69169475979376E-0311
8.69169475979376E-0311
8.69169475979376E-0311

但是对于 64 位，我得到了

 8.69169475979375E-0311
 0.00000000000000E+0000
 0.00000000000000E+0000

所以基本上 Tokyo 可以在 64 位模式下处理非正规数，常量被正确写入，但是从算术运算甚至使用 ldexp，非正规结果被刷新为零。

这个观察可以在其他系统上得到证实吗？如果是，它记录在哪里？ （我能找到的关于零刷新的唯一信息是， Denormals become zero when stored in a Real48)。

更新：我知道 32 位和 64 位都使用单个重载。对于 32 位，使用 x87 FPU，并且所有精度（单精度、双精度、扩展精度）的 ASM 代码几乎相同。 FPU 总是返回一个 80 位扩展，它存储在一个双精度数中，不会过早截断。 64 位代码在存储前进行精度调整。 同时我提交了一份问题报告 (https://quality.embarcadero.com/browse/RSP-20925)，重点关注 32 位或 64 位的不一致结果。

【参考方案1】：

更新：

只有编译器处理重载选择的方式有所不同。

正如@Graymatter 发现的那样，调用的LdExp 重载是32 位和64 位编译器的Single 类型。唯一的区别是代码库，其中 32 位编译器使用 asm 代码，而 64 位编译器使用 purepascal 实现。

要修复代码以使用正确的重载，请显式定义 LdExp() 第一个参数的类型，这样它可以工作（64 位）：

program Project116;

$APPTYPE CONSOLE
uses
  system.sysutils, system.math;
var
  x: double;
const
  twopm1030 : UInt64 = $0000100000000000; 2^(-1030)
begin
  x := PDouble(@twopm1030)^;
  writeln(x);
  x := ldexp(Double(1),-515);
  writeln(x*x);
  x := ldexp(Double(1),-1030);
  writeln(x);
  ReadLn;
end. 

输出：

 8.69169475979375E-0311
 8.69169475979375E-0311
 8.69169475979375E-0311

---

我会说这种行为应该报告为 RTL 错误，因为在您的案例中选择的重载函数是 Single 类型。结果类型是Double，编译器肯定应该做出相应的调整。 因为 32 位和 64 位编译器应该产生相同的结果。

---

注意，浮点类型的Double(1) 类型转换是在 Delphi 10.2 Tokyo 中引入的。有关先前版本中的解决方案，请参阅What is first version of Delphi which allows typecasts like double(10)。

【讨论】：

谢谢。我自己发现 ldexp 是在x := ldexp(1,-515) 之后打印 x 时的罪魁祸首。与 32 位相比，这至少是不一致的，并且独立于 $EXCESSPRECISION ON/OFF。我将提交问题报告。

这不是编译器错误。编译器在决定使用哪个重载时不会考虑返回值。这是 Ldexp 实现在使用 Win32 和 Win64 编译时给出不同结果的问题，因为它们的实现不同。 Win32 使用汇编程序，Win64 使用 pascal 代码。

不是编译器错误。

@Graymatter，谢谢。在 32 位编译器中，返回值在浮点寄存器中，它保留了 double 的精度，而 64 位编译器仅在 XMM 寄存器上以单精度运行。

这是 rtl 而不是编译器

【参考方案2】：

这里的问题是Ldexp(single) 会根据是否调用 ASM 代码或是否调用 pascal 代码返回不同的结果。在这两种情况下，编译器都会调用 Single 版本的重载，因为调用中没有指定类型。

您在 Win64 场景中执行的 pascal 代码尝试处理小于 -126 的指数，但该方法仍然无法正确计算结果，因为单个数字被限制为 8 位指数。汇编器似乎解决了这个问题，但我没有详细研究为什么会这样。

function Ldexp(const X: Single; const P: Integer): Single;
   Result := X * (2^P) 
$IFNDEF X86ASM
var
  T: Single;
  I: Integer;
const
  MaxExp =  127;
  MinExp = -126;
  FractionOfOne = $00800000;
begin
  T := X;
  Result := X;
  case T.SpecialType of
    fsDenormal,
    fsNDenormal,
    fsPositive,
    fsNegative:
      begin
        FClearExcept;
        I := P;
        if I > MaxExp then 
        begin
          T.BuildUp(False, FractionOfOne, MaxExp);
          Result := Result * T;
          I := I - MaxExp;
          if I > MaxExp then I := MaxExp;
        end
        else if I < MinExp then
        begin
          T.BuildUp(False, FractionOfOne, MinExp);
          Result := Result * T;
          I := I - MinExp;
          if I < MinExp then I := MinExp;
        end;
        if I <> 0 then
        begin
          T.BuildUp(False, FractionOfOne, I);
          Result := Result * T;
        end;
        FCheckExcept;
      end;
//    fsZero,
//    fsNZero,
//    fsInf,
//    fsNInf,
//    fsNaN:
  else
    ;
  end;
end;
$ELSE X86ASM
$IF     defined(CPUX86) and defined(ios) // iOS/Simulator
...
$ELSE 
asm // StackAlignSafe
        PUSH    EAX
        FILD    dword ptr [ESP]
        FLD     X
        FSCALE
        POP     EAX
        FSTP    ST(1)
        FWAIT
end;
$ENDIF
$ENDIF X86ASM

正如 LU RD 所建议的，您可以通过强制方法调用 Double 重载来解决此问题。有一个错误，但该错误是 ASM 代码与 Ldexp(const X: Single; const P: Integer) 中的 pascal 代码不匹配，而不是调用了不同的重载。