#include <boost/math/special_functions/gamma.hpp>
namespace boost{ namespace math{
template <class T>
calculated-result-type tgamma(T z);
template <class T, class Policy>
calculated-result-type tgamma(T z, const Policy&);
template <class T>
calculated-result-type tgamma1pm1(T dz);
template <class T, class Policy>
calculated-result-type tgamma1pm1(T dz, const Policy&);
}}
template <class T>
calculated-result-type tgamma(T z);
template <class T, class Policy>
calculated-result-type tgamma(T z, const Policy&);
Возвращает «истинную гамма» (отсюда и название tgamma) значения z:
Конечный аргументПолитикаявляется необязательным и может использоваться для контроля поведения функции: как она обрабатывает ошибки, какой уровень точности использовать и т. д. См. документациюдля более подробной информации.
Эффективно существует две версииtgammaфункции внутри: полностью общая версия, которая является медленной, но достаточно точной, и гораздо более эффективное приближение, которое используется там, где число цифр в значении T соответствует определенному приближениюLanczos. На практике любой построенный в плавающей точке тип, с которым вы столкнетесь, имеет соответствующееприближение Ланчоса, определенное для него. Также возможно, учитывая достаточное машинное время, генерировать дальнейшее приближениеLanczosс использованием программы libs/math/tools/lanczos_generator.cpp.
Тип возврата этой функции вычисляется с помощью.правила расчета типа результата:удваивается, когда T является целым типом, а T иначе.
template <class T>
calculated-result-type tgamma1pm1(T dz);
template <class T, class Policy>
calculated-result-type tgamma1pm1(T dz, const Policy&);
Возвращаетсяtgammadz+1-1. Внутренне реализация не использует сложение и вычитание, подразумеваемые определением, приводя к точным результатам даже для очень маленькихдз. Однако реализация ограничена либо 35-значной точностью, либо точностью приближенияLanczos, связанного с типом T, в зависимости от того, что более точно.
Тип возврата этой функции вычисляется с помощью.правила расчета типа результата:удваивается, когда T является целым типом, а T иначе.
Конечный аргументПолитикаявляется необязательным и может использоваться для контроля поведения функции: как она обрабатывает ошибки, какой уровень точности использовать и т. д. См. документациюдля более подробной информации.
В следующей таблице показаны пиковые ошибки (в единицах эпсилона), обнаруженные на различных платформах с различными типами плавающих точек, а также сравнения с другими распространенными библиотеками. Если не указано иное, любой тип плавающей точки, который является более узким, чем показанный, будет иметьфактически нулевую ошибку.
Table 6.1. Error rates for tgamma
|
|
Microsoft Visual C++ версия 12.0
Win32
|
GNU C++ версия 5.1.0
Linux
Double
|
GNU C++ версия 5.1.0
Linux
длинный двойной
|
Солнечный компилятор версии 0x5130
Солнечный солярис
|
факториалы |
Max = 1,85ε (Mean = 0,491ε)
:Max = 3,17ε (Mean = 0,928ε)] |
Max = 0ε (Mean = 0ε)
GSL 1.16:Max = 3.95ε (Mean = 0.783ε)]
Rmath 3.0.2:Max = 314ε (Mean = 93.4ε)]
Cephes:Max = 3.19ε (Mean = 0.884ε) |
Макс = 1,96ε (Mean = 0,483ε)
:Макс = 1,66ε (Mean = 0,584ε)]
:Макс = 1,66ε (Mean = 0,584ε)] |
Макс = 172ε (Mean = 41ε)
:Макс = 0ε (Mean = 0ε)] |
около 0 |
Макс = 1.96ε (Средний = 0.684ε)
:Макс = 1ε (Средний = 0.405ε)] |
Макс = 0ε (Mean = 0ε)
GSL 1.16:Макс = 4.51ε (Mean = 1,92ε)]
Rmath 3.0.2:Макс = 1ε (Mean = 0,335ε)]
Цефес:Макс = 1ε (Mean = 0,548ε)] |
Макс = 2ε (Mean = 0.73ε)
:Макс = 1ε (Mean = 0.376ε)]
:Макс = 1ε (Mean = 0.376ε)] |
Макс = 2ε (Средний = 0.647ε)
:Макс = 0,5ε (Средний = 0,0791ε)] |
около 1 |
Макс = 2ε (Mean = 0.865ε)
:Макс = 1ε (Mean = 0.4ε)] |
Макс = 0ε (Средний = 0ε)
GSL 1.16:Макс = 4.41ε (Средний = 1,81ε)]
Макс = 1ε Макс = 1ε [Средний = 0,518ε]] |
Макс = 2ε (Mean = 0.85ε)
:Макс = 0.918ε (Mean = 0.203ε)]
:Макс = 0.918ε (Mean = 0.203ε)] |
Макс = 3.01ε (Средний = 1,06ε)
:Макс = 1ε (Средний = 0,175ε)] |
около 2 |
Max = 2ε (Mean = 0.995ε)
:Max = 0ε (Mean = 0ε)] |
Макс = 0ε (Mean = 0ε)
GSL 1.16:Макс = 7.95ε (Mean = 3.12ε)]
Rmath 3.0.2:Макс = 1ε (Mean = 0,191ε)]
Цефес:Макс = 1,09ε (Mean = 0,502ε) |
Макс = 2ε (Mean = 0.913ε)
:Макс = 0,558ε (Mean = 0.298ε)]
:Макс = 0,558ε (Mean = 0.298ε) |
Макс = 5.01ε (Средний = 1,89ε)
:Макс = 0ε (Средний = 0ε)] |
около -10 |
Max = 1,73ε (Mean = 0,729ε)
:Max = 0,866ε (Mean = 0,445ε)] |
Макс = 0ε (Mean = 0ε)
GSL 1.16:Макс = 2.6ε (Mean = 1.05ε)]
Rmath 3.0.2:Макс = 6.34e+05ε (Mean = 1.2e+05ε)]
Цефес:Макс = 2.6ε (Mean = 0.956ε) |
Max = 2.6ε (Mean = 0.985ε)
:Max = 2.26ε (Mean = 1.08ε)]
:Max = 2.26ε (Mean = 1.08ε) |
Макс = 1,75ε (Средний = 0,819ε)
:Макс = 0ε (Средний = 0ε)] |
около -55 |
Max = 1.8ε (Mean = 0.817ε)
:Max = 3.87e+004ε (Mean = 6.71e+003ε)] |
Макс = 0ε (Mean = 0ε)
GSL 1.16:Макс = 1.8ε (Mean = 0.782ε)]
Rmath 3.0.2:Макс = 6.36e+06ε (Mean = 1.13e+06ε)]
Цефес:Макс = 2.7ε (Mean = 0.988ε) |
Макс = 1.8ε (Mean = 0.847ε)
:Макс = 1,79ε (Mean = 0,75ε)]
:Макс = 1,79ε (Mean = 0,75ε) |
Макс = 98.5ε (Средний = 53.4ε)
:Макс = 0ε (Средний = 0ε)] |
Table 6.2. Error rates for tgamma1pm1
|
|
Microsoft Visual C++ версия 12.0
Win32
|
GNU C++ версия 5.1.0
Linux
Double
|
GNU C++ версия 5.1.0
Linux
длинный двойной
|
Солнечный компилятор версии 0x5130
Солнечный солярис
|
tgamma1pm1(dz) |
Макс = 0,982ε (Средний = 0,399ε) |
Макс = 0ε (Средний = 0ε) |
Макс = 1,12ε (Средний = 0,49ε) |
Макс = 3.97ε (Средний = 0.713ε) |
Гамма относительно проста в тестировании: факториалы и полуцелочисленные факториалы могут быть рассчитаны точно другими способами и сопоставлены с гамма-функцией. Кроме того, некоторые тесты точности в известных сложных областях были рассчитаны с высокой точностью с использованием общей версии этой функции.
Функцияtgamma1pm1протестирована против значений, вычисленных очень наивно с использованием формулыtgamma1+dz—1с приближением lanczos, точным примерно до 100 десятичных цифр.
Общая версия функцииtgammaреализована в приближении Стерлинга к lgamma для больших z:
После экспоненциации используется нисходящая рекурсия для малых значений z.
Для типов известной точности используется аппроксимацияLanczos, класс признаковboost::математика::lanczos::lanczos_traitsкарты типа T с соответствующим приближением.
Для z в диапазоне -20< z< 1 затем используется рекурсия для перехода на z >1 через:
Для очень маленьких z это помогает сохранить идентичность:
Для z< -20 формула отражения:
используется. Особое внимание следует уделить оценкеz * sin(π * z): для уменьшения z перед умножением на π используется специальная процедура, обеспечивающая, чтобы результат находился в диапазоне [0, π/2]. Без этого в этом регионе происходит чрезмерное количество ошибок (что уже достаточно сложно, так как скорость изменения вблизи отрицательного полюса составляет).исключительновысоко.
Наконец, если аргумент представляет собой небольшое целое число, то используется поиск таблицы факториала.
Функцияtgamma1pm1реализована с использованием рациональных приближений, разработанных JMв области-0,5<dz<2. Это те же приближения (и внутренние процедуры), которые используются дляlgamma. Результатом приближения являетсяlogtgammadz+1], который может подаваться вexpm1для получения желаемого результата. За пределами диапазона-0,5<dz<2тогда наивная формулаtgamma1pm1dz=tgamma+11может использоваться непосредственно.
