跳至主要內容

第17章 数值

南风2024年10月26日大约 13 分钟

第17章 数值

  • 简介
  • 数学函数
  • 数值算法
    • 并行数值算法
  • 复数
  • 随机数
  • 向量算术
  • 数值极限
  • 类型别名
  • 数学常数
  • 建议

17.1 简介

C++的设计初衷并非主要针对数值计算。然而，数值计算通常发生在其他工作场景中——如科学计算、数据库访问、网络、仪器控制、图形、模拟和金融分析——因此，C++成为了执行更大系统中计算部分的一个有吸引力的工具。此外，数值方法已经远远超越了简单的对浮点数向量的循环处理。在计算中需要更复杂的数据结构时，C++的优势就变得相关了。其最终结果是，C++被广泛用于科学、工程、金融和其他涉及复杂数值计算的领域。因此，支持此类计算的功能和技术应运而生。本章将描述标准库中支持数值计算的部分。

17.2 数学函数

在 <cmath> 中，我们可以找到 标准数学函数 ，例如 sqrt() 、 log() 和 sin() ，它们适用于 float 、 double 和 long double 类型的参数：

abs(x)绝对值
ceil(x)大于等于 x 的最小整数
floor(x)小于等于 x 的最大整数
sqrt(x)平方根； x 必须非负
sin(x)正弦
cos(x)余弦
tan(x)正切
acos(x)反余弦；结果是非负的
asin(x)反正弦；返回最接近 0 的结果
atan(x)反正切
sinh(x)双曲正弦
cosh(x)双曲余弦
tanh(x)双曲正切
exp(x)自然底数 e 的指数
exp2(x)底数 2 的指数
log(x)自然对数，底数 e ； x 必须为正
log2(x)自然对数，底数 2 ； x 必须为正
log10(x)底数 10 的对数； x 必须为正

这些函数的复数版本（§17.4）在 <complex> 头文件中。对于每个函数，返回类型与参数类型相同。

通过设置 <cerrno> 中的 errno 为 EDOM 来报告定义域错误，为 ERANGE 来报告值域错误。例如：

errno = 0; // 清除旧的错误状态  
double d = sqrt(-1);  
if (errno == EDOM)  
    cerr << "sqrt()未定义负参数\n";  
errno = 0; // 清除旧的错误状态  
double dd = pow(numeric_limits<double>::max(), 2);  
if (errno == ERANGE)  
    cerr << "pow()的结果过大，无法表示为double类型\n";

更多的数学函数可以在 <cmath> 和 <cstdlib> 中找到。所谓的 特殊数学函数 ，如 beta() 、 rieman_zeta() 和 sph_bessel() ，也在 <cmath> 中。

17.3 数值算法

在 <numeric> 中，我们发现了一小组通用的数值算法，如 accumulate() 。

x=accumulate(b,e,i)x 是 i 和区间 [b:e) 内元素的总和
x=accumulate(b,e,i,f)使用 f 代替 + 来累积
x=inner_product(b,e,b2,i)x 是区间 [b:e) 和 [b2:b2+(e-b)) 的内积，即 i 加上每个 p1 （在 [b:e) 内）和对应的 p2 （在 [b2:b2+(e-b)) 内）的乘积之和
x=inner_product(b,e,b2,i,f,f2)使用 f 和 f2 代替加法和乘法来计算内积
p=partial_sum(b,e,out)区间 [out:p) 中第 i 个元素是 [b:b+i] 内元素的总和
p=partial_sum(b,e,out,f)使用 f 代替加法来计算部分和
p=adjacent_difference(b,e,out)区间 [out:p) 中第 i 个元素是 *(b+i)-*(b+i-1) ** （对于 ** i>0 ** ），如果 ** e-b>0 ** ，则 ** *out ** 是 ** *b
p=adjacent_difference(b,e,out,f)使用 f 代替减法来计算相邻差
iota(b,e,v)对于 [b:e) 内的每个元素，赋值为 v 并递增 ++v ；因此序列变为 v, v+1, v+2, ...
x=gcd(n,m)x 是整数 n 和 m 的最大公约数
x=lcm(n,m)x 是整数 n 和 m 的最小公倍数
x=midpoint(n,m)x 是 n 和 m 之间的中点

这些算法通过将常见操作（如计算总和）应用于各种序列，从而实现了操作的泛化。它们还将应用于这些序列元素的操作作为一个参数。对于每种算法，除了通用版本外，还提供了一个应用该算法最常见操作符的版本。例如：

list<double> lst {1, 2, 3, 4, 5, 9999.99999};  
auto s = accumulate(lst.begin(), lst.end(), 0.0);  // 计算总和：10014.9999

这些算法适用于所有标准库序列，并且可以将操作作为参数提供（§17.3）。

17.3.1 并行数值算法

在 <numeric> 中，数值算法（§17.3）具有与顺序版本略有不同的并行版本。特别是，并行版本允许以未指定的顺序对元素进行操作。并行数值算法可以采用执行策略参数（§13.6）： seq 、 unseq 、 par 和 par_unseq 。

x=reduce(b,e,v)x=accumulate(b,e,v) , 除了计算顺序不同
x=reduce(b,e)x=reduce(b,e,V{}) , 其中 V 是 b 的值类型
x=reduce(pol,b,e,v)x=reduce(b,e,v) 与执行策略 pol 结合使用
x=reduce(pol,b,e)x=reduce(pol,b,e,V{} , 其中 V 是 b 的值类型
p=exclusive_scan(pol,b,e,out)根据 pol , p=partial_sum(b,e,out) , 排除第 i 个输入元素在第 i 个累加中
p=inclusive_scan(pol,b,e,out)根据 pol , p=partial_sum(b,e,out) , 包含第 i 个输入元素在第 i 个累加中
p=transform_reduce(pol,b,e,f,v)对 [b:e) 范围内的每个 x 应用 f(x) , 然后进行 reduce 操作
p=transform_exclusive_scan(pol,b,e,out,f,v)对 [b:e) 范围内的每个 x 应用 f(x) , 然后进行 exclusive_scan 操作
p=transform_inclusive_scan(pol,b,e,out,f,v)对 [b:e) 范围内的每个 x 应用 f(x) , 然后进行 inclusive_scan 操作

为了简化，我省略了这些算法的以操作为参数的版本，而只使用了 + 和 = 。除了 reduce() 之外，我还省略了具有默认策略（顺序）和默认值的版本。

就像在 <algorithm> （§13.6）中的并行算法一样，我们可以指定一个执行策略：

vector<double> v {1, 2, 3, 4, 5, 9999.99999};  
auto s = reduce(v.begin(), v.end());  // 使用double类型作为累加器计算总和  
vector<double> large;  
// ... 用许多值填充large ...  
auto s2 = reduce(std::execution::par_unseq, large.begin(), large.end());  
// 使用可用的并行性计算总和

执行策略 par 、 seq 、 unseq 和 par_unseq 隐藏在 <execution> 中的 std::execution 命名空间中。

请测量并验证使用并行或向量化算法是否值得。

17.4 复数

标准库支持一系列复数类型，这些类型遵循了§5.2.1中描述的 complex 类的路线。为了支持标量为单精度浮点数（ float ）、双精度浮点数（ double ）等的复数，标准库中的 complex 是一个模板：

template<typename Scalar>  
class complex {  
public:  
    complex(const Scalar& re = {}, const Scalar& im = {});  // 默认函数参数	§3.4.1
    // ...  
};

对于复数，支持通常的算术运算和最常见的数学函数。例如：

void f(complex<float> fl, complex<double> db)  
{  
    complex<long double> ld {fl + sqrt(db)};  
    db += fl * 3;  
    fl = pow(1 / fl, 2);  
    // ...  
}

<complex> （§17.2）中定义了 sqrt() 和 pow() （指数运算）等通常的数学函数。

17.5 随机数

随机数在许多场合都非常有用，如测试、游戏、模拟和安全领域。 应用领域的多样性体现在标准库 <random> 提供的各种随机数生成器中。一个随机数生成器由两部分组成：

  1. 一个产生随机或伪随机值序列的 引擎
  2. 一个将这些值映射到某个范围内数学分布的 分布器

分布器的例子包括 uniform_int_distribution （产生的所有整数都是等可能的）、 normal_distribution （“钟形曲线”）和 exponential_distribution （指数增长）；每种分布都针对某个指定的范围。例如：

using my_engine = default_random_engine;  // 引擎类型 
using my_distribution = uniform_int_distribution<>; // 分布类型  
 
my_engine eng {};  // 引擎的默认版本 
my_distribution dist {1,6};  // 映射到整数1到6的分布
auto die = [&](){ return dist(eng); };  // 创建一个生成器
 
int x = die();  // 掷骰子：x成为一个[1:6]之间的值

由于其对通用性和性能的严格要求，一位专家将标准库的随机数组件视为“每个随机数库都想成为成长后的样子”。然而，它几乎不能被视为“对新手友好”。 using 语句和 lambda 表达式使得所做的事情更加明显。

对于新手（无论其背景如何），面向随机数库的完全通用接口可能构成一个严重的障碍。一个简单的均匀随机数生成器通常就足以开始。 例如：

Rand_int rnd {1,10};  // 创建一个[1:10]之间的随机数生成器 
int x = rnd();  // x是一个[1:10]之间的数

那么，我们如何实现这一点呢？我们必须得到一些像 die() 那样，在一个名为 Rand_int 的类中将引擎与分布器结合起来的东西：

class Rand_int {  
public:  
    Rand_int(int low, int high) :dist{low,high} { }  
    int operator()() { return dist(re); }  // 抽取一个整数  
    void seed(int s) { re.seed(s); }  // 选择新的随机引擎种子  
private:  
    default_random_engine re;  
    uniform_int_distribution<> dist;  
};

这个定义仍然是“专家级别”的，但在C++课程的第一周中，新手就可以学会使用 Rand_int() 。例如：

#include "pch.h";
#include <iostream>;
#include <random>

using namespace std;

class Rand_int {
public:
	Rand_int(int low, int high) :dist{ low,high } { }
	int operator()() { return dist(re); }  // 抽取一个整数  
	void seed(int s) { re.seed(s); }  // 选择新的随机引擎种子  
private:
	default_random_engine re;
	uniform_int_distribution<> dist;
};

int main()
{
	constexpr int max = 9;
	Rand_int rnd{ 0,max };  // 创建一个均匀随机数生成器  

	vector<int> histogram(max + 1);  // 创建一个适当大小的向量  
	for (int i = 0; i != 200; ++i)
		++histogram[rnd()];  // 用数字[0:max]的频率填充直方图 
	for (int i = 0; i != histogram.size(); ++i) { // 写出条形图  
		cout << i << '\t';
		for (int j = 0; j != histogram[i]; ++j) cout << '+';
		cout << '\n';
	}
}

输出是一个均匀分布（具有合理的统计变化）：

0       ++++++++++++++++
1       +++++++++++++++++++++
2       ++++++++++++++++++
3       ++++++++++++++++++++
4       ++++++++++++++++++
5       +++++++++++++++++
6       ++++++++++++++++
7       +++++++++++++++++++++++++++
8       +++++++++++++++++++
9       ++++++++++++++++++++++++++++

C++没有标准的图形库，所以我使用“ASCII图形”。显然，有很多开源和商业的图形和GUI库可用于C++，但在这本书中，我只使用ISO标准设施。

为了获得重复或不同的值序列，我们需要对引擎进行播种，即为其内部状态赋予一个新值。例如：

Rand_int rnd {10,20};  
for (int i = 0; i<10; ++i) cout << rnd() << ' ';  // 16 13 20 19 14 17 10 16 15 14
cout << '\n';  
rnd.seed(999);  
for (int i = 0; i<10; ++i) cout << rnd() << ' ';  // 11 17 14 19 20 13 20 14 16 19
cout << '\n';  
rnd.seed(999);  
for (int i = 0; i<10; ++i) cout << rnd() << ' ';  // 11 17 14 19 20 13 20 14 16 19
cout << '\n';

重复的序列对于确定性调试很重要。当我们不希望重复时，使用不同的值进行播种就很重要。如果你需要真正的随机数，而不是生成的伪随机序列，请查看你的机器上 random_device 的实现方式。

17.6 向量算术

第12.2节中描述的向量被设计为用于存储值的通用机制，它具有灵活性，并适应容器、迭代器和算法的架构。然而，它并不支持数学向量运算。向向量中添加此类运算很容易，但其通用性和灵活性排除了通常被认为对于严肃的数值工作至关重要的优化。因此，标准库在 <valarray> 中提供了一个类似向量的模板，称为 valarray ，它通用性较低，但更易于针对数值计算进行优化：

template<typename T>  
class valarray {  
// ...  
};

valarray 支持常用的算术运算和最常见的数学函数。例如：

void f(valarray<double>& a1, valarray<double>& a2)  
{  
    valarray<double> a = a1*3.14+a2/a1;  // 数值数组运算符*，+，/和=  
    a2 += a1*3.14;  
    a = abs(a);  
    double d = a2[7];  
    // ...  
}

这些运算是向量运算，即它们应用于所涉及向量的每个元素。

除了算术运算外， valarray 还提供步长访问以帮助实现多维计算。

17.7 数值极限

在 <limits> 中，标准库提供了描述内置类型属性的类，例如浮点数的最大指数或整数中的字节数。例如，我们可以断言 char 是有符号的：

static_assert(numeric_limits<char>::is_signed,"unsigned characters!");  
static_assert(100000<numeric_limits<int>::max(),"small ints!");

第二个断言有效，因为 numeric_limits<int>::max() 是一个 constexpr 函数（§1.6）。 我们可以为自定义的用户定义类型定义 numeric_limits 。

17.8 类型别名

基本类型（如 int 和 long long ）的大小是由实现定义的；也就是说，它们在不同的C++实现中可能不同。如果我们需要对整数的大小进行具体说明，可以使用 <stdint> 中定义的别名，如 int32_t 和 uint_least64_t 。后者表示至少有64位的无符号整数。

奇怪的后缀 _t 是C语言时代遗留下来的，当时认为名称应该反映其是一个别名。 其他常见的别名，如 size_t （ sizeof 运算符返回的类型）和 ptrdiff_t （从一个指针减去另一个指针的结果的类型），可以在 <stddef> 中找到。

17.9 数学常数

进行数学计算时，我们需要常见的数学常数，如 e 、 pi 和 log2e 。标准库提供了这些以及更多。它们有两种形式：一个允许我们指定确切类型的模板（例如， pi_v<T> ），以及最常见用法的简短名称（例如， pi 意味着 pi_v<double> ）。例如：

void area(float r)  
{  
using namespace std::numbers;  // 这里保存了数学常数  
double d = pi*r*r;  
float f = pi_v<float>*r*r;  
}

在这种情况下，差异很小（我们必须打印大约16位精度才能看到它），但在实际的物理计算中，这样的差异很快就会变得显著。在其他领域，如图形和人工智能，常数的精度也很重要，在这些领域中，较小的值表示变得越来越重要。

在 <numbers> 中，我们找到了 e （欧拉数）、 log2e （ e 的log2）、 log10e （ e 的log10）、 pi 、 inv_pi （ 1/pi ）、 inv_sqrtpi （ 1/sqrt(pi) ）、 ln2 、 ln10 、 sqrt2 （ sqrt(2) ）、 sqrt3 （ sqrt(3) ）、 inv_sqrt3 （ 1/sqrt3 ）、 egamma （欧拉-马斯切罗尼常数）和 phi （黄金比例）。当然，我们希望有更多不同领域的数学常数。 这很容易做到，因为这样的常数是变量模板，具有针对 double （或对于某个领域最有用的任何类型）的特化：

template<typename T>  
constexpr T tau_v = 2*pi_v<T>;  
constexpr double tau = tau_v<double>;

17.10 建议

  1. 数值问题往往很微妙。如果你对数值问题的数学方面不是百分之百确定，要么寻求专家建议，要么进行实验，或者两者都做；§17.1。
  2. 不要尝试仅使用基础语言进行严肃的数值计算；使用库；§17.1。
  3. 在编写循环从序列中计算值之前，考虑使用 accumulate() 、 inner_product() 、 partial_sum() 和 adjacent_difference() ；§17.3。
  4. 对于大量数据，尝试使用并行和向量化算法；§17.3.1。
  5. 使用 std::complex 进行复数算术；§17.4。
  6. 将引擎绑定到分布以获得随机数生成器；§17.5。
  7. 注意你的随机数对于你的预期用途是否足够随机；§17.5。
  8. 不要使用C标准库的 rand() ；它对于实际用途来说不够随机；§17.5。
  9. 当运行时效率比操作和元素类型的灵活性更重要时，使用 valarray 进行数值计算；§17.6。
  10. 数值类型的属性可通过 numeric_limits 访问；§17.7。
  11. 使用 numeric_limits 检查数值类型是否适合其用途；§17.7。
  12. 如果你想明确整数类型的大小，请为其使用别名；§17.8。
上次编辑于: 2024/11/8 04:35:03
贡献者: nfzwy
默认页脚
Copyright © 2024 南风