В различных контекстах, например, для сокращения аргумента для математических функций, необходимо вычислить(a - K) / (a + K), гдеa является аргументом положительной переменной иK постоянная Во многих случаях,K является степенью двойки, которая является вариантом использования, относящимся к моей работе. Я ищу эффективные способы вычисления этого коэффициента более точно, чем это можно сделать с помощью простого деления. Можно предположить аппаратную поддержку слияния с многократным добавлением (FMA), поскольку в настоящее время эта операция предоставляется всеми основными архитектурами ЦП и ГП и доступна в C / C ++ через функцииfma() а такжеfmaf().

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

 /* Compute q = (a - K) / (a + K) with improved accuracy. Variant 1 */
 m = a - K;
 p = a + K;
 r = 1.0f / p;
 q = m * r;
 t = fmaf (q, -2.0f*K, m);
 e = fmaf (q, -m, t);
 q = fmaf (r, e, q);

Для аргументовa в промежутке[K/2, 4.23*K]Приведенный выше код вычисляет частное почти правильно округленное для всех входных данных (максимальная ошибка очень близка к 0,5 ульпс), при условии, чтоK является степенью 2, и в промежуточных результатах нет переполнения или недостатка. ЗаK этот код не является степенью двойки, но он все же более точен, чем простой алгоритм, основанный на делении. С точки зрения производительности этот код может бытьБыстрее чем наивный подход на платформах, где обратные числа с плавающей точкой могут быть вычислены быстрее, чем деление с плавающей точкой.

Я делаю следующее наблюдение, когдаK = 2n: Когда верхняя граница рабочего интервала увеличивается до8*K, 16*K, ... максимальная ошибка постепенно увеличивается и начинает медленно приближаться к максимальной ошибке наивного вычисления снизу. К сожалению, то же самое не кажется верным для нижней границы интервала. Если нижняя граница падает до0.25*Kмаксимальная погрешность улучшенного метода выше равна максимальной погрешности наивного метода.

Существует ли метод для вычисления q = (a - K) / (a ​​+ K), который может обеспечить меньшую максимальную ошибку (измеряется вULP по сравнению с математическим результатом) по сравнению как с наивным методом, так и с вышеуказанной кодовой последовательностью, в более широком интервале,в частности, для интервалов, нижняя граница которых меньше0.5*K? Эффективность важна, но скорее всего допустимо несколько больше операций, чем используется в приведенном выше коде.

В одном ответе ниже было указано, что я мог бы повысить точность, возвращая частное в виде неоцененной суммы двух операндов, то есть в виде пары голова-хвостq:qloто есть похож на хорошо известныйfloat и дваждыdouble форматы. В моем коде выше это означало бы изменение последней строки наqlo = r * e.

Этот подход, безусловно, полезен, и я уже обдумывал его использование для логарифма с повышенной точностью для использования вpow(), Но это принципиально не помогает с желаемым расширением интервала, в котором расширенные вычисления обеспечивают более точные коэффициенты. В конкретном случае я смотрю, я хотел бы использоватьK=2 (для одинарной точности) илиK=4 (для двойной точности), чтобы сохранить интервал первичного приближения узким, а интервал дляa примерно [0,28]. Практическая проблема, с которой я сталкиваюсь, состоит в том, что для аргументов <0,25 * K точность улучшенного деления не намного лучше, чем при использовании наивного метода.