高性能蒙特卡洛积分模拟框架
4.31/5 (11投票s)
用于金融衍生品定价的高性能蒙特卡洛模拟框架。
引言
在本文中，我编写了一个高性能的蒙特卡洛模拟框架。为了实现高性能，我同时使用了并行进程和每个进程中的并行多线程。MPI（消息传递接口）协议用于协调多个进程（本文使用 openmpi 实现的 MPI），C++ 11 线程库用于实现多线程。由于我来自金融行业，这个蒙特卡洛框架被设计用于对金融衍生品进行定价。我使用静态继承，特别是“好奇的递归模板模式”（curiously recurring template pattern），使框架足够灵活，可以为具有不同类型支付的任何类型的衍生品定价，同时保持高性能，并避免在运行时动态继承时可能产生的虚函数开销。
多进程和多线程
蒙特卡洛模拟非常适合并行计算，因为所有路径都相互独立，并且我们需要运行大量路径才能收敛。为了实现高性能，我使用了两个层次的并行。第一层是并行进程。MPI 接口提供了进程集之间的通信。首先，我将模拟路径的计算分成不同的进程，然后我继续在每个进程中将路径分成多个线程。在所有线程完成（使用线程 join 来等待所有线程）后，我将这些结果组合起来，作为该进程的最终结果。在所有进程完成后，我使用 MPI reduce 操作将来自进程的结果组合起来，以获得蒙特卡洛模拟的最终结果。
MPI 代码和蒙特卡洛模拟的初始化如下
int main(int argc, char** argv)
{
std::chrono::duration<double> elapsed_seconds_0;
auto start_0 = std::chrono::high_resolution_clock::now();
/* MPI Init*/
int proc_id = 0;
int num_procs = 1;
int namelen;
float result_each_proc;
float result_sum_global;
char proc_name[MPI_MAX_PROCESSOR_NAME];
MPI_Init( &argc, &argv );
MPI_Comm_size( MPI_COMM_WORLD, &num_procs );
MPI_Comm_rank( MPI_COMM_WORLD, &proc_id );
MPI_Get_processor_name( pro_name, &namelen );
int num_path_per_process = SIM_NUM_PATH / num_procs;
Instrument Stock(100, 0.5, 0.01);
//EuropeanCallOptionPayOff CallPayOff(100);
//DerivativeInstrument<EuropeanCallOptionPayOff> DerivedInstr(0.5, &Stock, &CallPayOff);
//MonteCarlo<EuropeanCallOptionPayOff> MonteCalc(STEP_NUM, num_path_per_process, &DerivedInstr,NUM_THREAD);
// EuropeanPutOptionPayOff PutPayOff(100);
// DerivativeInstrument<EuropeanPutOptionPayOff> DerivedInstr(0.5, &Stock, &PutPayOff);
// MonteCarlo<EuropeanPutOptionPayOff> MonteCalc(STEP_NUM, num_path_per_process, &DerivedInstr,NUM_THREAD);
// EuropeanBarrierDownOutCallOptionPayOff BarrierOutPayOff( 100, 100);
// DerivativeInstrument<EuropeanBarrierDownOutCallOptionPayOff> DerivedInstr(0.5, &Stock, &BarrierOutPayOff);
// MonteCarlo<EuropeanBarrierDownOutCallOptionPayOff> MonteCalc(STEP_NUM, num_path_per_process, &DerivedInstr,NUM_THREAD);
//EuropeanBarrierDownInCallOptionPayOff BarrierInPayOff( 100, 100);
//DerivativeInstrument<EuropeanBarrierDownInCallOptionPayOff> DerivedInstr(0.5, &Stock, &BarrierInPayOff);
//MonteCarlo<EuropeanBarrierDownInCallOptionPayOff> MonteCalc(STEP_NUM, num_path_per_process, &DerivedInstr,NUM_THREAD);
AsianCallOptionPayOff AsianPayOff( 100 );
DerivativeInstrument<AsianCallOptionPayOff> DerivedInstr(0.5, &Stock, &AsianPayOff);
MonteCarlo<AsianCallOptionPayOff> MonteCalc(STEP_NUM, num_path_per_process, &DerivedInstr,NUM_THREAD);
std::chrono::duration<double> elapsed_seconds;
auto start = std::chrono::high_resolution_clock::now();
elapsed_seconds_0 = start - start_0;
double dRes = MonteCalc.CalcMultiThreadPlainMC();
MPI_Allreduce ( &result_each_proc, &result_sum_global, 1, MPI_FLOAT, MPI_SUM, MPI_COMM_WORLd );
float final_result = result_sum_global / num_procs;
auto end = std::chrono::high_resolution_clock::now();
elapsed_seconds = end - start;
if ( 0 == proc_id )
{
std::cout << "MPI Initial" << num_procs << " processes spend" << elapsed_seconds_0.count() << "second" << std::endl;
std::cout << SIM_NUM_PATH << "number of " << STEP_NUM << "steps of simulation takes " << elapsed_seconds.count()<<"second" << std::endl;
}
MPI_Finalize();
return 0;
下表显示了在一台 2.2GHz 的 32 核机器上，使用不同数量的进程和线程所花费的时间。
Thread/Process 1 2 5 10 20 32
1 149.176 74.785 29.726 14.916 10.351 10.136
2 74.698 37.351 15.561 10.419 8.538 7.591
5 29.803 15.040 9.500 7.609 7.530 6.926
10 11.967 9.035 7.298 7.008 6.850 6.757
16 12.470 9.454 7.151 6.822 6.880 6.758
从性能表中，我们可以得出结论，增加进程数量或线程数量都可以提高性能，直到达到某个限制。主机是 32 核机器，我们无法仅通过使用 32 个进程就达到瓶颈，并且增加线程肯定会提高性能。此外，MPI 多进程的初始化时间会随着进程数量的增加而增加。
下表显示了在一台 2.2GHz 的 32 核机器上，使用不同数量的进程和线程所花费的时间。
Process 1 2 5 10 20 32
Time to
initialize processes 0.049 0.052 0.063 0.165 0.855 1.973
蒙特卡洛模拟框架
蒙特卡洛模拟是一类广泛的计算算法，它依赖于重复的随机抽样来获得数值结果；通常，我们会进行大量的模拟以获得未知概率实体的分布。蒙特卡洛模拟可用于数值积分以计算期望平均值，而本框架旨在计算金融衍生品价格的期望平均值。本文中的代码示例使用一个结构体来表示随机标的资产价格作为状态变量。该框架的用户当然可以在此结构体中添加更多变量，以扩展到具有多个随机过程和特定相关方差的多维度。我使用 State 的 vector 来定义一条路径，并使用 paths 的 vector 来定义来自蒙特卡洛模拟的所有路径。
struct State
{
State( float price ) : m_price ( price ) {}
float m_price;
};
typedef std::vector< State > Path;
此处编码的随机标的资产过程遵循具有无风险利率漂移项和随机项的几何布朗运动。随机项规定，与前一状态价格的百分比变化来自均值为零、方差为标的资产价格方差乘以时间变化的正态分布的随机样本。
for (int num = 0; num < lSimNum; num++)
{
pVecPath = &(m_vecPaths[pathindex+num]);
(*pVecPath)[0].m_price = InitialPrice;
S = InitialPrice;
fPathSum = 0;
for ( int i = 1; i < m_iSteps; i++ )
{
bPathValid = true;
number = distribution(generator);
volterm = voldt* number;
S = S * exp(drift + volterm);
(*pVecPath)[i].m_price = (*pVecPath)[i-1].m_price * exp( drift + volterm );
}
}
不同种类的金融衍生品具有不同类型的支付。我在这里编码的支付类型包括与路径无关的期权（如香草看涨期权或看跌期权）、路径相关的障碍期权（如障碍期权），以及沿路径的平均标的资产价格（如亚式期权）。我使用静态继承，特别是“好奇的递归模板模式”，使框架足够灵活，可以为具有不同类型支付的任何类型的衍生品定价，同时保持高性能，并避免在运行时动态继承时可能产生的虚函数开销。
基础支付类（base payoff class）的实际实现如下
template < class Derived>
class PayOff
{
public:
float CalcPayOff(float dUnderPrice)
{
return static_cast <Derived *> (this)->CalcPayOff(dUnderPrice);
}
bool CheckPathBarrier( float fPrice )
{
return static_cast < Derived * >(this)->CheckPathBarrier( fPrice );
}
void CalcSumPath( float fPrice, float * fRet )
{
return static_cast< Derived * > (this)->CalcSumPath( fPrice, fRet );
}
PayOffType m_type;
};
以下代码显示了欧洲看涨期权如何继承支付基类。
class EuropeanCallOptionPayOff : public PayOff <EuropeanCallOptionPayOff >
{
public:
EuropeanCallOptionPayOff(float Strike) : m_strike(Strike)
{
m_type = PathIndependent;
}
float CalcPayOff(float dPrice)
{
if (dPrice > m_strike)
return dPrice - m_strike;
else
return 0;
}
private:
float m_strike;
};
Derivative 类是一个模板类，以 payoff 作为模板参数。代码如下
template<typename DerivatedPayOff>
class DerivativeInstrument
{
public:
DerivativeInstrument(float Time, Instrument * pIns, PayOff<DerivatedPayOff> * pPO) : dTimeToMature(Time), pInstrument(pIns), pPayOff(pPO)
{
}
float dTimeToMature;
Instrument * pInstrument;
PayOff<DerivatedPayOff> * pPayOff;
};
关注点
下载代码时，您可以下载 openmpi 的头文件和库，并在 makefile 中正确设置它们，以便代码正常工作。如果您有任何问题，请告诉我。
历史
2014/08/31，初次发布
2014/09/01，将代码添加到本文中
2014/10/12，修改论文和代码以泛化此框架