本征ConditionType数组:一种有效的广播方式,而不是循环

问题描述

我有一段对性能至关重要的代码,我需要检查一个数组中低于阈值的值,然后有条件地设置其他两个数组的值。我的代码如下:

#include <Eigen/Dense>

int main(){
    Eigen::ArrayXXd
        a (1,100),b (2,c (3,100);
    
    a.setRandom();
    b.setRandom();
    c.setRandom();
    
    constexpr double minVal { 1e-8 };
    
    /* the code segment in question */
    /* option 1 */
    for ( int i=0; i<2; ++i ){
        b.row(i)   = (a < minVal).select( 0,c.row(i+1) / a );
        c.row(i+1) = (a < minVal).select( 0,c.row(i+1) );
    }
    /* option 2,which is slower */
    b = (a < minVal).replicate(2,1).select( 0,c.bottomRows(2) / a.replicate(2,1) );
    c.bottomRows(2) = (a < minVal).replicate(2,c.bottomRows(2) );

    return 0;
}

检查其值是否达到阈值a的数组minVal具有一行并且具有动态的列数。其他两个数组bc分别具有两行和三行,并且列数与a相同。

现在,我想以一种更eigen的方式执行上述逻辑,而在选项1中没有该循环,因为通常eigen会提高性能,我从不希望在编写原始循环时匹配。 但是,我唯一想到的方法是选项2,它明显比选项1慢。

做上述事情的正确而有效的方法是什么?还是循环已经是我最好的选择?

解决方法

您可以尝试以下操作:

  • 用固定的行数和动态的列数定义数组类型,即,您可以将 Eigen :: ArrayXXd 替换为 Eigen :: Array
  • 使用固定大小的块操作(请参见https://eigen.tuxfamily.org/dox/group__TutorialBlockOperations.html),即,您可以将 bottomRows(N)替换为 bottomRows ()并类似地 replicate(2,1) replicate ()

我已经更改了代码中的数组类型,并包括了我提到的可能的改进的第三个选项:

#include <Eigen/Dense>

#include <iostream>
#include <chrono>

constexpr int numberOfTrials = 1000000;
constexpr double minVal{ 1e-8 };

typedef Eigen::Array<double,1,Eigen::Dynamic> Array1Xd;
typedef Eigen::Array<double,2,Eigen::Dynamic> Array2Xd;
typedef Eigen::Array<double,3,Eigen::Dynamic> Array3Xd;

inline void option1(const Array1Xd& a,Array2Xd& b,Array3Xd& c)
{
    for (int i = 0; i < 2; ++i) {
        b.row(i) = (a < minVal).select(0,c.row(i + 1) / a);
        c.row(i + 1) = (a < minVal).select(0,c.row(i + 1));
    }
}

inline void option2(const Array1Xd& a,Array3Xd& c)
{
    b = (a < minVal).replicate(2,1).select(0,c.bottomRows(2) / a.replicate(2,1));
    c.bottomRows(2) = (a < minVal).replicate(2,c.bottomRows(2));
}

inline void option3(const Array1Xd& a,Array3Xd& c)
{
    b = (a < minVal).replicate<2,1>().select(0,c.bottomRows<2>() / a.replicate<2,1>());
    c.bottomRows<2>() = (a < minVal).replicate<2,c.bottomRows<2>());
}

int main() {
    Array1Xd a(1,100);
    Array2Xd b(2,100);
    Array3Xd c(3,100);

    a.setRandom();
    b.setRandom();
    c.setRandom();

    auto tpBegin1 = std::chrono::steady_clock::now();
    for (int i = 0; i < numberOfTrials; i++)
        option1(a,b,c);
    auto tpEnd1 = std::chrono::steady_clock::now();

    auto tpBegin2 = std::chrono::steady_clock::now();
    for (int i = 0; i < numberOfTrials; i++)
        option2(a,c);
    auto tpEnd2 = std::chrono::steady_clock::now();

    auto tpBegin3 = std::chrono::steady_clock::now();
    for (int i = 0; i < numberOfTrials; i++)
        option3(a,c);
    auto tpEnd3 = std::chrono::steady_clock::now();

    std::cout << "(Option 1) Average execution time: " << std::chrono::duration_cast<std::chrono::microseconds>(tpEnd1 - tpBegin1).count() / (long double)(numberOfTrials) << " us" << std::endl;
    std::cout << "(Option 2) Average execution time: " << std::chrono::duration_cast<std::chrono::microseconds>(tpEnd2 - tpBegin2).count() / (long double)(numberOfTrials) << " us" << std::endl;
    std::cout << "(Option 3) Average execution time: " << std::chrono::duration_cast<std::chrono::microseconds>(tpEnd3 - tpBegin3).count() / (long double)(numberOfTrials) << " us" << std::endl;

    return 0;
}

我获得的平均执行时间如下(i7-9700K,msvc2019,启用优化,NDEBUG):

(Option 1) Average execution time: 0.527717 us
(Option 2) Average execution time: 3.25618 us
(Option 3) Average execution time: 0.512029 us

并启用了AVX2 + OpenMP:

(Option 1) Average execution time: 0.374309 us
(Option 2) Average execution time: 3.31356 us
(Option 3) Average execution time: 0.260551 us

我不确定这是否是最“本征”的方法,但我希望它能有所帮助!

相关问答

依赖报错 idea导入项目后依赖报错,解决方案:https://blog....
错误1:代码生成器依赖和mybatis依赖冲突 启动项目时报错如下...
错误1:gradle项目控制台输出为乱码 # 解决方案:https://bl...
错误还原:在查询的过程中,传入的workType为0时,该条件不起...
报错如下,gcc版本太低 ^ server.c:5346:31: 错误:‘struct...