近期关注本个人号的朋友们可能要失望了，虽然本人工作日几乎每天都要“算命”，但公众号中却很少写算命的内容。一是，由于合规等原因，不方便讲市场情况，更不能讲个券，二是，折腾这个公众号的原因主要是周末在泡图书馆的时候找点放松的事情做。
春节的时候，启动了一个研究的题目，但由于节日气氛浓厚的缘故，做成了烂尾，近期才继续。本着如果不写一点东西，就完全没有动力的缘故，后面择机、逐期将一些结果、心得等写下来。还是老规矩，如果没有涉及到具体品种等容易有利益冲突的内容，就写详细一些；如果涉及到容易与工作利益冲突的内容，就“打码”。前面几篇，由于内容都很简单，所有代码就都不“打码”了。
一、欧式期权
研究任务是算一些比较复杂的奇异期权的价值，选用的方法是Monte Carlo Simulation。它具有很容易上手的优良特点，当然，计算速度比较慢，所以比较适合于研究Python高性能计算。
在研究复杂的期权结构之前，先从欧式开始，原因是它有解析解，可以用来验证方法的正确性。
Black-Scholes-Merton公式的推导就不讲了，直接上Python代码：
![]()
前面一段是的import各类库，对python不太熟，这段花了我不少时间找函数对应的库名称，可能找到的不是速度最优的，但能用就行。用Matlab例子所用的参数，调用这个函数，最后可以计算得到call和put价值。
![]()
得到的结果是，Call=13.695272738608146， Put=6.349714381299741，与Matlab官网的结果基本一致，说明这个函数没有写错。
后面那个N，是蒙特卡洛模拟时用的参数，这里暂时不用的。
二、蒙特卡洛模拟
欧式期权只在期末才是行权日，所以实际上只要模拟期末的标的资产价格就行。但这里我们选择把资产价格的整条路径都模拟出来，主要是为了便于后期计算路径依赖型衍生品价格。
模拟很简单，就是设定N条路径，每一条路径从期初股价S0开始，T+1日的股价是T日股价的一点点扰动，扰动服从正态分布，分布的均值是 r-σ*σ/2, 标准差为σ。这里先用最直观的循环代码写出来。
![]()
这段代码非常直观，而且非常省内存空间（每条路径每天股价只需要一个变量储存、每条路径的结果加总也只储存一个变量）。但如果代码写成这样，那么就完全丧失了让程序自动优化的可能了，因为里面每一个值都是生成下一个值得前提条件（同一条路径上前后两期的值、相邻两条路径的值等）。
这个函数有多慢？我们测试可知。测试用的机器是Surface Pro 6， CPU是i5 8250U，性能应该还算勉强够用的。但即使这样，我也不想计算N=300万条路径，而是选择N/10，即30万条。调用代码是：
- (call0, put0)=mcopt_serial(S0,K,r,T,sig,N//10)
模拟的结果与精确值都有一些误差，这个误差会随着模拟路径增加而有所下降。具体误差计算放在后文。
PS：如上所述，这段代码写得实在太“实诚”，我测试用numba.jit 去做加速，完全没有效果。也许是参数不太对，如果有哪位读者有jit加速好的参数欢迎回复后台。
三、挖掘CPU速度潜能
这里掉书袋，深入浅出讲讲CPU的性能。CPU是电脑的计算能力源泉，常听到说2.4GHz、3.6GHz、4.5GHz这样的数字，描述的是CPU的频率，简单可以理解为一秒钟计算多少次，例如1GHz就是一秒计算10亿次，2GHz就是20亿次每秒。这个数据当然是越高就越快。但遗憾的是，电流速度是极限，现在常用的CPU到5GHz已经差不多顶天了。
除了加快频率，CPU还从多个方面发展。
打个比方，你可以将CPU看做大学生，为了让它计算更快，一种方法就是如上所说的加快频率，即要这个大学生提高速度，但这种方法很容易就到头，频率越高，功率上升特别快，很容易达到系统极限；
第二种方法是，多找几个大学生一起来算，这就是所谓的多核心了。最初CPU里面只有单核，后面发现加一个核心，可能计算速度就达到2倍，随着CPU制程升级，核心可以做得更小，于是就可以加更多的核心。现在（1903），X东自营上能买到的最多核心的CPU是线程撕裂者2990WX，它有32个核心。但这种方法的弊端是，如果一件事情只能由一个人做，那么多余的那些核心可能就只能在旁边干瞪眼，所以它只适合可以多个人一起协作的工作。
第三种方法，让核心更高效率的工作。想象这个大学生计算完一个任务后，需要去拿数据进行下一次计算，如果系统能够在当前任务还在进行之际，就将下一个任务交到这个人手上，那么速度就能比原来有所提高。这种方法叫做超线程技术。相当于让一个核心一次性地先后做两个任务。例如，前述的2990WX，每个核心有2个线程，即一共64线程。这种方法的弊端就是，虽然看上去是一个核心当做两个用，但性能增长受到任务性质的影响，对于大多数任务而言，提速效果很不明显。
第四种方法，对核心计算能力进行升级，里面又分为两类，一是优化计算的流程，相当于从大一升级到大四，做事情效率能提高；二是增加特殊的计算模块，更快速执行一些计算任务。
后者更值得一说。早期多媒体计算需求刚刚兴起，但CPU计算能力弱，于是Intel在CPU中加了一个模块，专门计算多媒体数据，调用这个模块需要用到特殊的指令，这些指令就被称为MMX指令集。
后来，大家发现这样做，成本低，效果好。于是，一发不可收拾，又搞了SSE家族、AVX家族等指令集，主要目标是一次性多计算几个数据。
例如，我现在电脑8250U支持AVX2指令集，它一次性可以计算256位/4个双精度浮点数据，性能比CPU常规计算流程快很多。可以用CPU-Z等工具查看你的CPU性能数据。
![]()
四、Python中充分利用CPU的高性能
现在CPU几乎都是多核，但我们常用的Cpython由于GIL的缘故，自己搞并行计算很麻烦。但没有太大的必要，如果你用的Python版本里面自带了MKL库（Anaconda近年的版本中自带了），那么只要你程序多用Numpy等系统自带的科学计算包，它就能比较充分地利用CPU多核心、多种类指令集的优点。
这个里面的调用逻辑是： 你的程序——>Numpy/Scipy等库——>MKL——>CPU
当然，Python用Numpy可以加速还有一大原因是，Python作为一种解释性语言，循环等语句本身效率就很低，改为系统自带的包后相当于直接底层调用。
五、Numpy版本的程序
如上所述，我们要提高速度，就要充分利用Numpy自带的功能，当然同时研究减少了自己写的循环代码，即所谓代码的向量化。通过查找相关函数，我们可以将前面代码中的两段循环用一个语句计算出来。
![]()
其中paths就是计算得到的N条股价路径，每一条有T日的数据。对于欧式期权而言，不需要这么复杂，但之后计算路径依赖型期权时就用得上了。此外，这段代码还计算了结果的标准差，它可以用来衡量模拟所得结果的精确度。一般而言当模拟路径增长到原来得4倍时，精确度可以提升到原来的2倍。
调用代码如下：
- (call1, put1, std1)=mcopt(S0,K,r,T,sig,N)
耗费时间是9.34秒，与前面代码计算300万条预计耗费560秒相比，速度提高了很多很多。而且，如果观察计算过程中CPU的占用，也可以发现，多个核心都没有闲着。例如我现在这台电脑是4核心8线程，CPU占用50%算是用满所有核心。
![]()
PS：这段代码有一个较大的问题，很耗内存。如果要模拟更多的路径，就需要换一台更多内存的电脑，或者分批进行了。
六、Intel和AMD的CPU之争
AMD的CPU价格便宜量又足（核心多），可是做科学计算的时候却比较少用它（超算中心不管，主要考虑家用级），那么为什么呢？
我个人找到两个原因，第一原因是截止到当前，AMD的CPU对AVX等科学计算敏感的支持不够，有人研究发现AMD的CPU，直到Zen架构，它的AVX指令的Vector Width/Units，都显著不及同期Intel的CPU，自然比不过。今年下半年发布的Zen2架构倒是显著改进了，希望到时候能有性能大幅提升。
![]()
（图来源：http://numberworld.org/y-cruncher/faq.html）
第二个原因就是MKL本身。大家看看MKL的全称就知道了：Intel Math Kernel Library，开发人员明确说了，他们只在Intel的CPU上做了测试。
本篇闲扯到此结束，下一篇准备写CUDA，即显卡计算。我们常规的计算都是在CPU上进行的，而下一篇的计算是要搬到显卡上进行，就是如下图这样的硬件：
|
|
转播
