C++ 中的快速百分位数

我的程序计算风险价值指标的蒙特卡罗模拟。为了尽可能简化，我有：

1/ simulated daily cashflows
2/ to get a sample of a possible 1-year cashflow, 
   I need to draw 365 random daily cashflows and sum them

因此，每日现金流量是根据经验给出的分布函数，需要采样 365 次。为此，我

 1/ sort the daily cashflows into an array called *this->distro*
 2/ calculate 365 percentiles corresponding to random probabilities

我需要对年度现金流进行这种模拟，例如 10K 次，才能获得模拟年度现金流的总体数据。准备好每日现金流量的分配函数后，我进行抽样...

for ( unsigned int idxSim = 0; idxSim < _g.xSimulationCount; idxSim++ )
{
    generatedVal = 0.0;
    for ( register unsigned int idxDay = 0; idxDay < 365; idxDay ++ )
    {
        prob = (FLT_TYPE)fastrand();         // prob [0,1]
        dIdx = prob * dMaxDistroIndex;       // scale prob to distro function size
                                             // to get an index into distro array
        _floor = ((FLT_TYPE)(long)dIdx);     // fast version of floor
        _ceil  = _floor + 1.0f;              // 'fast' ceil:)
        iIdx1  = (unsigned int)( _floor );
        iIdx2  = iIdx1 + 1;

        // interpolation per se
        generatedVal += this->distro[iIdx1]*(_ceil - dIdx  );
        generatedVal += this->distro[iIdx2]*(dIdx  - _floor);
    }
    this->yearlyCashflows[idxSim] = generatedVal ;
}

两者里面的代码forCycles 进行线性插值。如果说 1000 美元对应 prob=0.01，10000 美元对应 prob=0.1，那么如果我没有 p=0.05 的经验数，我想通过插值得到 5000 美元。

问题是：这段代码运行正确，尽管探查器表示该程序将 cca 60% 的运行时间花在了插值本身上。所以我的问题是，如何才能更快地完成这项任务？ VTune 报告的特定线路的示例运行时间如下：

prob = (FLT_TYPE)fastrand();         //  0.727s
dIdx = prob * dMaxDistroIndex;       //  1.435s
_floor = ((FLT_TYPE)(long)dIdx);     //  0.718s
_ceil  = _floor + 1.0f;              //    -

iIdx1  = (unsigned int)( _floor );   // 4.949s
iIdx2  = iIdx1 + 1;                  //    -

// interpolation per se
generatedVal += this->distro[iIdx1]*(_ceil - dIdx  );  //    -
generatedVal += this->distro[iIdx2]*(dIdx  - _floor);  // 12.704s

破折号表示探查器报告这些行没有运行时间。

任何提示将不胜感激。 丹尼尔

EDIT:c.fogelklou 和 MSalters 都指出了很大的改进。符合 c.fogelklou 所说的最佳代码是

converter = distroDimension / (FLT_TYPE)(RAND_MAX + 1)
for ( unsigned int idxSim = 0; idxSim < _g.xSimulationCount; idxSim++ )
{
    generatedVal = 0.0;
    for ( register unsigned int idxDay = 0; idxDay < 365; idxDay ++ )
    {
        dIdx  = (FLT_TYPE)fastrand() * converter;
        iIdx1 = (unsigned long)dIdx);
        _floor = (FLT_TYPE)iIdx1;
        generatedVal += this->distro[iIdx1] + this->diffs[iIdx1] *(dIdx  - _floor);
    }
}

我所拥有的最好的 MSalter 的路线是

normalizer = 1.0/(FLT_TYPE)(RAND_MAX + 1);
for ( unsigned int idxSim = 0; idxSim < _g.xSimulationCount; idxSim++ )
{
    generatedVal = 0.0;
    for ( register unsigned int idxDay = 0; idxDay < 365; idxDay ++ )
    {
        dIdx  = (FLT_TYPE)fastrand()* normalizer ;
        iIdx1 = fastrand() % _g.xDayCount;
        generatedVal += this->distro[iIdx1];
        generatedVal += this->diffs[iIdx1]*dIdx;
    }
}

第二个代码大约是。速度提高 30%。现在，在 95 秒的总运行时间中，最后一行消耗了 68 秒。最后一行仅消耗 3.2 秒，因此 double*double 乘法一定是魔鬼。我想到了 SSE - 将最后三个操作数保存到一个数组中，然后执行 this->diffs[i]*dIdx[i] 的向量乘法并将其添加到 this->distro[i] 但此代码运行了 50%慢点。因此，我想我碰壁了。

非常感谢大家。 D .

这是一项小型优化的提议，消除了对 ceil、两次强制转换和一次乘法的需要。如果您在定点处理器上运行，这可以解释为什么 float 和 int 之间的乘法和转换花费如此长的时间。在这种情况下，请尝试使用定点优化或在编译器中打开浮点（如果 CPU 支持）！

for ( unsigned int idxSim = 0; idxSim < _g.xSimulationCount; idxSim++ )
{
    generatedVal = 0.0;
    for ( register unsigned int idxDay = 0; idxDay < 365; idxDay ++ )
    {
        prob = (FLT_TYPE)fastrand();         // prob [0,1]
        dIdx = prob * dMaxDistroIndex;       // scale prob to distro function size
                                             // to get an index into distro array
        iIdx1  = (long)dIdx;
        _floor = (FLT_TYPE)iIdx1;     // fast version of floor
        iIdx2  = iIdx1 + 1;

        // interpolation per se
        {
           const FLT_TYPE diff = this->distro[iIdx2] - this->distro[iIdx1];
           const FLT_TYPE interp = this->distro[iIdx1] + diff * (dIdx - _floor);
           generatedVal += interp;
        }
    }
    this->yearlyCashflows[idxSim] = generatedVal ;
}