Boost C++ 库

伯努利数 是一个有理数序列，可用于泰勒级数展开、欧拉-麦克劳林公式和黎曼 zeta 函数。

伯努利数用于评估一些 Boost.Math 函数，包括 tgamma、lgamma 和 polygamma 函数。

单个伯努利数

概要

#include <boost/math/special_functions/bernoulli.hpp>

namespace boost { namespace math {

template <class T>
T bernoulli_b2n(const int n);  // Single Bernoulli number (default policy).

template <class T, class Policy>
T bernoulli_b2n(const int n, const Policy &pol); // User policy for errors etc.

}} // namespaces

描述

两者都返回第 (2 * n)^th 伯努利数 B_2n。

请注意，由于所有奇数编号的伯努利数都为零（除了 B₁ 是 -½），因此该接口将仅返回偶数编号的伯努利数。

此函数对低索引伯努利数使用快速表查找，而较大的值则根据需要计算并缓存。缓存机制需要一定量的线程安全代码，因此 unchecked_bernoulli_b2n 可能为性能关键代码提供更好的接口。

最后的 Policy 参数是可选的，可用于控制函数的行为：它如何处理错误、要使用的精度级别等。

有关更多详细信息，请参阅 Policies。

示例

一个简单的示例计算 B₄ 的值，其中返回类型为 double，请注意 bernoulli_b2n 的参数是 2 而不是 4，因为它计算 B_2N。

try
{ // It is always wise to use try'n'catch blocks around Boost.Math functions
  // so that any informative error messages can be displayed in the catch block.
std::cout
  << std::setprecision(std::numeric_limits<double>::digits10)
  << boost::math::bernoulli_b2n<double>(2) << std::endl;

因此 B₄ == -1/30 == -0.0333333333333333

如果我们使用 Boost.Multiprecision 及其 50 位十进制数字浮点类型 cpp_dec_float_50，我们可以计算更大的数字（如 B₂₀₀）的值，并获得更高的精度。

std::cout
  << std::setprecision(std::numeric_limits<boost::multiprecision::cpp_dec_float_50>::digits10)
  << boost::math::bernoulli_b2n<boost::multiprecision::cpp_dec_float_50>(100) << std::endl;

-3.6470772645191354362138308865549944904868234686191e+215

单个（未检查）伯努利数

概要

#include <boost/math/special_functions/bernoulli.hpp>

template <>
struct max_bernoulli_b2n<T>;

template<class T>
inline T unchecked_bernoulli_b2n(unsigned n);

unchecked_bernoulli_b2n 提供对伯努利数的访问，不进行任何溢出或无效参数检查。它被实现为直接（且非常快速）的表查找，虽然不建议用于一般用途，但在需要最终性能且错误检查移到循环外部的内部循环中可能很有用。

您可以传递给 ... 的最大值是 ...：传递大于该值的值将导致缓冲区溢出错误，因此在使用此直接接口时，在您自己的代码中放置错误处理显然很重要。

... 的值因类型 T 而异，对于类型 ...，它是不会溢出目标类型的最大值：例如，... 是 129。但是，对于多精度类型，它是结果可以表示为两个 64 位整数之比的最大值，例如 ... 仅为 17。当然，更大的索引可以传递给 ...，但那样您将失去快速表查找（即值可能需要计算）。

/*For example:
*/
   std::cout << "boost::math::max_bernoulli_b2n<float>::value = "  << boost::math::max_bernoulli_b2n<float>::value << std::endl;
   std::cout << "Maximum Bernoulli number using float is " << boost::math::bernoulli_b2n<float>( boost::math::max_bernoulli_b2n<float>::value) << std::endl;
   std::cout << "boost::math::max_bernoulli_b2n<double>::value = "  << boost::math::max_bernoulli_b2n<double>::value << std::endl;
   std::cout << "Maximum Bernoulli number using double is " << boost::math::bernoulli_b2n<double>( boost::math::max_bernoulli_b2n<double>::value) << std::endl;

boost::math::max_bernoulli_b2n<float>::value = 32
Maximum Bernoulli number using float is -2.0938e+038
boost::math::max_bernoulli_b2n<double>::value = 129
Maximum Bernoulli number using double is 1.33528e+306

多个伯努利数

概要

#include <boost/math/special_functions/bernoulli.hpp>

namespace boost { namespace math {

// Multiple Bernoulli numbers (default policy).
template <class T, class OutputIterator>
OutputIterator bernoulli_b2n(
  int start_index,
  unsigned number_of_bernoullis_b2n,
  OutputIterator out_it);

// Multiple Bernoulli numbers (user policy).
template <class T, class OutputIterator, class Policy>
OutputIterator bernoulli_b2n(
  int start_index,
  unsigned number_of_bernoullis_b2n,
  OutputIterator out_it,
  const Policy& pol);
}} // namespaces

描述

提供了两个版本的伯努利数函数，用于在一个调用中计算多个伯努利数（一个使用默认策略，另一个允许用户定义的策略）。

这些返回一系列伯努利数

示例

我们可以通过一次调用计算并保存所有浮点精度的伯努利数。

std::vector<float> bn; // Space for 32-bit `float` precision Bernoulli numbers.

// Start with Bernoulli number 0.
boost::math::bernoulli_b2n<float>(0, 32, std::back_inserter(bn)); // Fill vector with even Bernoulli numbers.

for(size_t i = 0; i < bn.size(); i++)
{ // Show vector of even Bernoulli numbers, showing all significant decimal digits.
    std::cout << std::setprecision(std::numeric_limits<float>::digits10)
        << i*2 << ' '
        << bn[i]
        << std::endl;
}

0 1
2 0.166667
4 -0.0333333
6 0.0238095
8 -0.0333333
10 0.0757576
12 -0.253114
14 1.16667
16 -7.09216
18 54.9712
20 -529.124
22 6192.12
24 -86580.3
26 1.42552e+006
28 -2.72982e+007
30 6.01581e+008
32 -1.51163e+010
34 4.29615e+011
36 -1.37117e+013
38 4.88332e+014
40 -1.92966e+016
42 8.41693e+017
44 -4.03381e+019
46 2.11507e+021
48 -1.20866e+023
50 7.50087e+024
52 -5.03878e+026
54 3.65288e+028
56 -2.84988e+030
58 2.38654e+032
60 -2.14e+034
62 2.0501e+036

当然，对于任何浮点类型，都存在一个最大伯努利数，该数在溢出指数之前可以计算出来。默认策略下，如果我们尝试计算过高的伯努利数，则会抛出异常。

try
{
  std::cout
  << std::setprecision(std::numeric_limits<float>::digits10)
  << "Bernoulli number " << 33 * 2 <<std::endl;

  std::cout << boost::math::bernoulli_b2n<float>(33) << std::endl;
}
catch (std::exception ex)
{
  std::cout << "Thrown Exception caught: " << ex.what() << std::endl;
}

我们将收到有用的错误消息（前提是使用了 try'n'catch 块）。

Bernoulli number 66
Thrown Exception caught: Error in function boost::math::bernoulli_b2n<float>(n):
Overflow evaluating function at 33

此示例的源代码位于 bernoulli_example.cpp

精度

所有函数通常返回的值都在浮点类型的一个 ULP（最后一位单位）内。

实现

实现细节在 bernoulli_details.hpp 和 unchecked_bernoulli.hpp 中。

对于 ...，这是通过从静态初始化的表进行简单表查找来实现的；对于较大的 ... 值，这是通过正切数算法实现的，如论文中所述：《伯努利数、正切数和正割数的快速计算》，Richard P. Brent 和 David Harvey，http://arxiv.org/pdf/1108.0286v3.pdf (2011)。

正切（或 Zag）数（偶数交替排列数）被定义，并且其中也给出了它们的生成函数。

正切数与伯努利数 B_i 的关系由 Brent 和 Harvey 的公式 14 给出

它们与伯努利数 B_i 的关系定义为

如果 i > 0 且 i 为偶数，则

否则如果 i == 0，则 B_i = 1
否则如果 i == 1，则 B_i = -1/2
否则如果 i < 0 或 i 为奇数，则 B_i = 0

请注意，计算值存储在固定大小的表中，通过原子操作（即无锁编程）进行线程安全访问，这使得对缓存值的访问开销远低于预期 - 通常对于多精度类型，线程同步的成本可以忽略不计，而对于内置类型，此代码通常无论如何都不会执行。对于无法合理计算或存储在我们缓存中的非常大的参数，使用了 Luschny 的渐近展开式