Использование 80-битного FPU в Delphi и Pascal для повышения точности вычислений с плавающей запятой

Delphi , Синтаксис , Память и Указатели

В мире программирования, где точность вычислений играет ключевую роль, особенно в научных и финансовых приложениях, возникает необходимость в более высокой точности представления чисел с плавающей запятой. Стандартный тип Double (или Real в Pascal) предоставляет 64-битную точность, что может быть недостаточно для некоторых задач. В этой статье мы рассмотрим возможности использования 80-битного FPU (Floating-Point Unit) в Delphi и Pascal для повышения точности вычислений.

Проблема: Ограниченность 64-битной точности

64-битное представление чисел с плавающей запятой (Double) имеет свои ограничения. Оно может привести к ошибкам округления и потере точности при выполнении сложных вычислений, особенно когда речь идет о большом количестве итераций или работе с очень большими или очень маленькими числами. Это может существенно повлиять на результаты моделирования, финансовых расчетов и других критически важных приложений.

Решение: Использование 80-битного FPU

80-битное FPU предоставляет значительно большую точность, чем 64-битное. Это позволяет уменьшить ошибки округления и получить более точные результаты вычислений. Однако, использование 80-битного FPU не всегда простое и требует определенных усилий.

Как это сделать в Delphi и Pascal?

Недавно, обсуждение на форуме Lazarus/Free Pascal (https://forum.lazarus.freepascal.org/index.php/topic,62198.0.html) выявило, что использование 80-битного FPU в Pascal может быть осуществлено с помощью:

1. Использование sfpux80 и ufloatx80: Первоначально предлагалось использовать эти юниты. Однако, они не всегда доступны "из коробки" и могут потребовать компиляции из исходного кода. По сообщениям пользователей, эти юниты доступны в исходном коде FPC (Free Pascal Compiler) и могут быть найдены в rtl/inc.

2. Использование кода Stephen L. Moshier: Более практичным и доступным решением является использование публичного кода Stephen L. Moshier (http://www.moshier.net/#Rubid_pc). Он предоставляет библиотеки для работы с 80-битными числами с плавающей запятой. Этот код требует компиляции с использованием MinGW (Minimalist GNU for Windows). Для компиляции требуется включить заголовки <string.h> и <malloc.h>.

Пример кода на Object Pascal (Delphi) с использованием кода Stephen L. Moshier:

{$MODE OBJFPC}{$H+}
{$MINFPCONSTPREC 64}
{$LINK ioldoubl}
{$LINKLIB libmsvcrt}
{$LINKLIB libgcc}

program ioldouble_test1;
uses ctypes, sysutils, windows, math;

type
  extfloat = record
    ext_number : array[0..15] of uint8;
  end;
type
  extfloatp = ^extfloat;

function strtold ( s : pchar; var ps : pchar):extfloat; cdecl; external name '_strtold';
procedure ldtostr(x:extfloatp; string_: pchar; ndigs:uint32; flags:int32; fmt: pchar);cdecl; external name '_IO_ldtostr';
function ldtoa(d:extfloatp; mode:int32; ndigits:int32; var decpt:int32; var sign:int32; var rve:pchar):pchar; cdecl; external name '_IO_ldtoa';
function snprintf ( s:pchar; n:csize_t; format:pchar):cint; varargs; cdecl; external name 'snprintf';
function atof(s: PChar):double; cdecl; external name 'atof';

function extsqrt(x: extfloat):extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fldt     x
      fsqrt
      fstpt y
    end;
    extsqrt:=y;
end;

function extadd(a: extfloat; b: extfloat):extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fldt     a
      fldt     b
      faddp   %st, %st(1)
      fstpt y
    end;
    extadd:=y;
end;

function extsub(a: extfloat; b: extfloat):extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fldt     a
      fldt     b
      fsubrp   %st, %st(1)
      fstpt y
    end;
    extsub:=y;
end;

function extmul(a: extfloat; b: extfloat):extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fldt     a
      fldt     b
      fmulp   %st, %st(1)
      fstpt y
    end;
    extmul:=y;
end;

function extdiv(a: extfloat; b: extfloat):extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fldt     a
      fldt     b
      fdivrp   %st, %st(1)
      fstpt y
    end;
    extdiv:=y;
end;

function extval( s:ansistring):extfloat;
var
  pc, ppc:pchar;
  y:extfloat;
begin
  pc:=pchar(s+#0);
  ppc:=#0;
  y:=strtold(pc, ppc);
  extval:=y;
end;

function extstr( x:extfloat):ansistring;
var
  s:pchar;
  res:ansistring;
begin
  res:='                                                                        ';
  s:=pchar(res+#0);
  ldtostr(@x, s, 17, 0, 'E');
  extstr:=s;
end;

function extstr( x:double):ansistring;
var
  s:pchar;
  res:ansistring;
begin
  res:='                                                                        ';
  s:=pchar(res+#0);
  snprintf(s, 64, '%.15g', x);
  extstr:=s;
end;

operator := (r : ansistring) z : double;
var s : pchar;
begin
  s:=pchar(r);
  z:=atof(s);
end;

operator := (r : ansistring) z : extfloat;
begin
  z:=extval(r);
end;

operator := (r : int32) z : extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fild    r
      fstpt y
    end;
  z:=y;
end;

operator := (r : int64) z : extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fildl    r
      fstpt y
    end;
  z:=y;
end;

operator := (r : double) z : extfloat;
var
  y:extfloat;
begin
  {$ASMMODE att}
    asm
      fldl    r
      fstpt y
    end;
  z:=y;
end;

operator + (x : extfloat; y : extfloat) z : extfloat;
begin
  z:=extadd(x, y);
end;

operator - (x : extfloat; y : extfloat) z : extfloat;
begin
  z:=extsub(x, y);
end;

operator * (x : extfloat; y : extfloat) z : extfloat;
begin
  z:=extmul(x, y);
end;

operator / (x : extfloat; y : extfloat) z : extfloat;
begin
  z:=extdiv(x, y);
end;

function sqrt(x : extfloat): extfloat;Overload;
begin
  result:=extsqrt(x);
end;

// main test code
var x, y, z:extfloat;
begin
  x:='1.23456789';
  y:=123;
  z:=x+y;
  writeln(extstr(z));
  x:=8;
  y:=9;
  z:=x/y;
  writeln(extstr(z));
  x:=2;
  z:=sqrt(x);
  writeln(extstr(z));
end.

Альтернативные решения:

• Использование библиотеки libmingwex: Эта библиотека предоставляет готовые функции для работы с 80-битным FPU, что упрощает разработку. Она также требует компиляции с использованием MinGW.
• Использование библиотеки GMP (GNU Multiple Precision Arithmetic Library): GMP - это высокопроизводительная библиотека для работы с числами произвольной точности. Хотя она не использует FPU напрямую, она позволяет достичь очень высокой точности вычислений. Интеграция GMP с Delphi/Pascal может быть сложной, но предоставляет максимальную гибкость.

Заключение:

Использование 80-битного FPU в Delphi и Pascal позволяет значительно повысить точность вычислений, что особенно важно для приложений, требующих высокой точности. Хотя процесс может быть немного сложным, особенно при использовании кода Stephen L. Moshier или библиотеки libmingwex, результат стоит затраченных усилий. Выбор конкретного метода зависит от ваших потребностей и уровня сложности, который вы готовы принять. Помните, что правильный выбор типа данных и алгоритмов может значительно повысить производительность и точность ваших программ.

Создано по материалам из источника по ссылке.

Статья рассматривает методы повышения точности вычислений в Delphi и Pascal за счёт использования 80-битного FPU, включая примеры кода и альтернативные решения.

Комментарии и вопросы

Получайте свежие новости и обновления по Object Pascal, Delphi и Lazarus прямо в свой смартфон. Подпишитесь на наш Telegram-канал delphi_kansoftware и будьте в курсе последних тенденций в разработке под Linux, Windows, Android и iOS

:: Главная :: Память и Указатели ::