数值算法与随机数：从枯燥计算到蒙特卡洛模拟

数值算法这块，说实话，很多C++开发者平时用得不多。但一旦遇到性能敏感的数据处理任务，你会发现STL提供的这些工具简直是救命稻草。我个人习惯在写数据分析类代码时，优先考虑用<numeric>里的算法，而不是自己手写循环。

今天咱们就聊聊数值算法和随机数生成。嗯，这两个话题其实紧密相关——你想想看，蒙特卡洛模拟就是它们结合的典型场景。

accumulate：最基础的累加器

std::accumulate可能是<numeric>里最常用的函数了。它的作用很简单：从一个初始值开始，依次对序列中的每个元素执行二元操作。

#include <numeric>
#include <vector>
#include <iostream>

int main() {
    std::vector<int> data = {1, 2, 3, 4, 5};
    
    // 默认加法：1+2+3+4+5 = 15
    int sum = std::accumulate(data.begin(), data.end(), 0);
    std::cout << "Sum: " << sum << std::endl;
    
    // 自定义操作：乘法
    int product = std::accumulate(data.begin(), data.end(), 1, 
                                  std::multiplies<int>());
    std::cout << "Product: " << product << std::endl;
    
    return 0;
}

inner_product：内积运算

说白了，std::inner_product就是两个序列对应元素相乘再求和。这在向量运算、相似度计算里非常常见。

#include <numeric>
#include <vector>

std::vector<double> a = {1.0, 2.0, 3.0};
std::vector<double> b = {4.0, 5.0, 6.0};

// 1*4 + 2*5 + 3*6 = 32
double dot = std::inner_product(a.begin(), a.end(), b.begin(), 0.0);

你还可以提供自定义的加法和乘法操作，实现更复杂的组合逻辑。比如计算两个字符串序列的某种匹配得分。

partial_sum 与 adjacent_difference

这两个函数是一对「好兄弟」。一个做前缀和，一个做差分。

#include <numeric>
#include <vector>
#include <iterator>

std::vector<int> v = {1, 2, 3, 4, 5};
std::vector<int> result(5);

// 前缀和：1, 3, 6, 10, 15
std::partial_sum(v.begin(), v.end(), result.begin());

// 差分：1, 1, 1, 1, 1
std::adjacent_difference(v.begin(), v.end(), result.begin());

iota：连续序列生成器

std::iota的名字有点怪，但功能很简单：从起始值开始，每次递增1，填充整个范围。

#include <numeric>
#include <vector>

std::vector<int> indices(10);
std::iota(indices.begin(), indices.end(), 0);
// indices = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9}

嗯，这里要注意：iota不会检查溢出。如果你给的类型是uint8_t，从250开始，很快就会出现回绕。我在嵌入式项目中吃过这个亏。

随机数引擎与分布

C++11之后，随机数生成有了全新的设计。核心思想是：引擎负责产生随机比特，分布负责把这些比特转换成你需要的分布形式。

随机数引擎

常用的引擎有：

• std::mt19937：梅森旋转算法，质量好、速度快，我90%的场景都用它
• std::default_random_engine：实现定义，不推荐跨平台使用
• std::random_device：真随机数（如果硬件支持），适合用来做种子

#include <random>

// 推荐做法：用random_device生成种子
std::random_device rd;
std::mt19937 gen(rd());

分布

分布决定了随机数的统计特性：

分布	用途
uniform_int_distribution	整数均匀分布，比如掷骰子
uniform_real_distribution	实数均匀分布，比如[0,1)随机数
normal_distribution	正态分布，模拟测量误差
bernoulli_distribution	伯努利分布，抛硬币

std::uniform_real_distribution<double> dist(0.0, 1.0);
double r = dist(gen);  // 生成[0,1)之间的随机数

实战：蒙特卡洛模拟估算π

好了，理论说完了，咱们来点实战。蒙特卡洛模拟的核心思想是：用大量随机采样来逼近某个数值结果。

估算π的经典方法：在一个正方形内随机撒点，统计落在内切圆中的比例。这个比例乘以4就是π的近似值。

#include <iostream>
#include <random>
#include <cmath>

double estimate_pi(int num_samples) {
    std::random_device rd;
    std::mt19937 gen(rd());
    std::uniform_real_distribution<double> dist(-1.0, 1.0);
    
    int inside = 0;
    for (int i = 0; i < num_samples; ++i) {
        double x = dist(gen);
        double y = dist(gen);
        if (x*x + y*y <= 1.0) {
            ++inside;
        }
    }
    
    return 4.0 * inside / num_samples;
}

int main() {
    std::cout << "π ≈ " << estimate_pi(1000000) << std::endl;
    return 0;
}

运行结果大概是3.1415左右。采样越多，精度越高。这就是蒙特卡洛的魅力——用随机性解决确定性问题。

知识体系总览

下面这张图帮你理清本章的核心脉络：

数值算法和随机数，看似是两个独立的话题，但在蒙特卡洛模拟中完美融合。你掌握了这些工具，就能高效地处理大量数值计算任务。

记住一点：STL提供的这些算法，不仅仅是帮你少写几行代码。它们经过了严格的性能优化和正确性验证。我在项目中用partial_sum替代手写循环后，性能提升了将近30%。这就是标准库的价值。

---

公众号:蓝海资料掘金营，微信deep3321