上回把三重循环换了下顺序，矩阵乘法从9秒干到1.5秒，当时对照的OpenBLAS成绩是0.05秒，中间还隔着三十倍，这三十倍总不能全赖人家多线程吧。翻了一堆资料把这笔账算明白了，大头其实藏在两个地方：CPU的指令级并行，和矩阵的分解策略。
写代码的人很容易有个惯性思维：一行 a[i]*b[i] 就是做了一次乘法，CPU老老实实算完一个再算下一个。这是指令串行的想象，真实情况是现代CPU一条SIMD指令一次就能算一批：AVX2的寄存器256位宽，塞得下4个double，一条乘加指令（FMA）出发就是4乘4加；再加上主流CPU有两个FMA运算端口，一个时钟周期理论吞吐是16个浮点操作。标量代码再怎么优化，每个周期也就吃掉这一小口，理论峰值的大餐一口没动。
有些问题Deterministic的路子根本走不通，硬解是天文数字，但是拿随机数去"扔"，反而能又快又好的逼近答案。这就是蒙特卡洛方法，高性能计算领域里最气质独特的一支。
最经典的入门例子，往单位正方形里随机撒点，数一数落在四分之一圆里的比例，命中数 $m$ 与总样本数 $N$ 之比就是面积的估计，于是 $\hat{\pi} = \dfrac{4m}{N}$：
起因是同事说同样三重循环算矩阵乘法，只是把循环变量的顺序换了一下，快了几十倍，一开始我是不信的，算法都没换哪来的几十倍。自己跑了一遍被啪啪打脸，顺手把背后的原理捋清楚了，故有此文备忘。
C语言里教科书式的矩阵乘法，N阶方阵要做 $2N^3$ 次浮点运算，一行注释都不用多加：
在研究如何利用thrift开发高扩展、高并发、高容错的rpc框架的学习过程中有了本篇文章做备忘。例子中使用的通信方式在进行进一步封装扩展后就完全可以用于生产环境。本例并非常见的简单例子。