как выполнять высокоточные вычисления в D?

Для некоторых университетских работ мне нужно аппроксимировать некоторые числа - например, число Эйлера с рядом. Поэтому я должен добавить очень маленькие числа, но у меня проблемы с точностью. Если число очень маленькое, это не влияет на результат.

real s;  //sum of all previous terms
ulong k; //factorial

s += 1.0/ k;

после каждого шага k становится еще меньше, но после 10-го раунда результат уже не меняется и застрял на 2,71828.


person NaN    schedule 09.02.2011    source источник


Ответы (3)


Если вам нужно решение, которое будет работать с использованием собственных типов, вы сможете получить разумные результаты, пытаясь всегда добавлять числа одинаковой величины. Один из способов сделать это — вычислить первые X членов ряда, а затем несколько раз заменить два наименьших числа этой суммой:

auto data = real[N];
foreach(i, ref v; data) {
  v = Fn(i);
}

while(data.length > 1) {
  data.sort(); // IIRC .sort is deprecated but I forget what replaced it.
  data[1] += data[0];
  data = data[1..$];
}

return data[0];

(Мин куча сделает это немного быстрее.)

person BCS    schedule 10.02.2011
comment
отлично.. попробую так! - person NaN; 10.02.2011

Типы с плавающей запятой с фиксированной точностью, которые изначально поддерживаются модулем с плавающей запятой вашего процессора (float, double, real), не оптимальны для любых вычислений, требующих многозначной точности, таких как приведенный вами пример.

Проблема в том, что эти типы с плавающей запятой имеют конечное число цифр точности (фактически двоичных цифр), что ограничивает длину числа, которое может быть представлено таким типом данных. Тип float имеет ограничение примерно в 7 десятичных цифр (например, 3.141593); тип double ограничен 14 (например, 3.1415926535898); и тип real имеет аналогичный лимит (чуть больше, чем у double).

Таким образом, добавление чрезвычайно малых чисел к значению с плавающей запятой приведет к потере этих цифр. Посмотрите, что произойдет, если мы сложим вместе следующие два значения с плавающей запятой:

float a = 1.234567f, b = 0.0000000001234567
float c = a + b;

writefln("a = %f b = %f c = %f", a, b, c);

И a, и b являются действительными значениями с плавающей запятой и сохраняют примерно 7 цифр точности каждый по отдельности. Но при добавлении сохраняются только самые первые 7 цифр, потому что они возвращаются обратно в число с плавающей запятой:

1.2345670001234567 => 1.234567|0001234567 => 1.234567
                              ^^^^^^^^^^^
                         sent to the bit bucket

Таким образом, c оказывается равным a, потому что более мелкие цифры точности от сложения a и b отбрасываются.

Вот еще одно объяснение концепции, вероятно, намного лучше моего.


Ответом на эту проблему является арифметика произвольной точности. К сожалению, поддержка арифметики с произвольной точностью не реализована в оборудовании ЦП; поэтому его (обычно) нет на вашем языке программирования. Однако существует множество библиотек, поддерживающих типы с плавающей запятой произвольной точности и необходимые для них математические операции. См. этот вопрос для некоторых предложений. Вероятно, сегодня вы не найдете каких-либо библиотек, специфичных для D, для этой цели, но существует множество библиотек C (GMP, MPFR и т. д.), которые должны быть достаточно просты в использовании по отдельности, и тем более, если вы сможете найти D привязки для одного из них.

person jgottula    schedule 09.02.2011

Как уже упоминалось, вам нужно использовать какую-то стороннюю библиотеку арифметики с плавающей запятой с высокой точностью (я думаю, что у Tango или Phobos есть только модуль для целочисленной арифметики произвольной длины).

dil — это проект D, использующий MPFR. Там вы должны найти привязки.

person Trass3r    schedule 10.02.2011