C++ 基于元素的代数表达式的基本表达式模板

示例


介绍和动机


表达式模板(以下称为ET)是一种强大的模板元编程技术,用于加快有时非常昂贵的表达式的计算速度。它广泛用于不同领域,例如在线性代数库的实现中。

对于此示例,请考虑线性代数计算的上下文。更具体地说,仅涉及按元素运算的计算。这种计算是ET的最基本应用,它们很好地介绍了ET内部的工作方式。

让我们看一个激励性的例子。考虑表达式的计算:

Vector vec_1, vec_2, vec_3;

// 初始化vec_1,vec_2和vec_3。

Vector result = vec_1 + vec_2*vec_3;

在这里为了简单起见,我假设类Vector和运算符+(向量加:逐个元素的加法运算)和运算符*(这里表示向量内积:也逐个元素的运算)都正确实现,在数学上应该是这样。

在不使用ET(或其他类似技术)的常规实现中,为了获得最终结果,至少进行了五种Vector实例的构造result:

  1. 对应于三个实例vec_1,vec_2和vec_3。

  2. 代表的结果的临时Vector实例。_tmp_tmp = vec_2*vec_3;

  3. 最后,通过正确使用返回值优化,最终构建resultin result = vec_1 + _tmp;。

使用ET的实现可以消除 Vector _tmp在2中创建临时对象,从而仅留下实例的四种构造Vector。更有趣的是,考虑以下更复杂的表达式:

Vector result = vec_1 + (vec_2*vec3 + vec_1)*(vec_2 + vec_3*vec_1);

Vector实例总共将有四种构造:vec_1, vec_2, vec_3和result。换句话说,在此示例中,仅涉及基于元素的操作,因此可以确保不会从中间计算中创建任何临时对象


ET如何工作


基本上,用于任何代数计算的ET包含两个构建块:

  1. 纯代数表达式PAE):它们是代数表达式的代理/抽象。纯代数不执行实际的计算,它们只是计算工作流的抽象/建模。PAE可以是任何代数表达式的输入或输出的模型。PAE实例通常被认为是廉价复制的。

  2. 惰性评估:这是实际计算的实现。在下面的示例中,我们将看到对于仅涉及元素操作的表达式,惰性求值可以在最终结果的索引访问操作内部实现实际计算,从而创建了按需求值的方案:不执行计算仅当访问/要求最终结果时。

那么,在此示例中,我们具体如何实现ET?让我们现在来看一下。

始终考虑以下代码片段:

Vector vec_1, vec_2, vec_3;

// 初始化vec_1,vec_2和vec_3。

Vector result = vec_1 + vec_2*vec_3;

计算结果的表达式可以进一步分解为两个子表达式:

  1. 向量加表达式(表示为plus_expr

  2. 向量内积表达式(表示为innerprod_expr)。

什么外星人做的是以下情况:

  • ET不会立即计算每个子表达式,而是先使用图形结构为整个表达式建模。图中的每个节点代表一个PAE。节点的边缘连接代表实际的计算流程。因此,对于上面的表达式,我们获得下图:

           result = plus_expr( vec_1, innerprod_expr(vec_2, vec_3) )
             /   \
            /     \
           /       \
          /   innerprod_expr( vec_2, vec_3 )
         /         /  \
        /         /    \
       /         /      \
    vec_1     vec_2    vec_3
  • 最终的计算是通过查看图的层次结构实现的:由于这里我们仅处理基于元素的运算,因此每个索引值in的计算result 都可以独立完成:的最终评估result可以延迟地推迟到基于元素的评估的每个元素result。换句话说,由于元素的计算result,elem_res可以使用对应于元素表示vec_1(elem_1), (vec_2)elem_2和vec_3(elem_3)为:

    elem_res = elem_1 + elem_2*elem_3;

因此,无需创建临时Vector存储中间内部乘积的结果:一个元素的整个计算可以完全一起完成,并在索引访问操作中进行编码


这是运行中的示例代码。


File vec.hh:std :: vector的包装器,用于在调用构造时显示日志。


#ifndef EXPR_VEC
# define EXPR_VEC

# include <vector>
# include <cassert>
# include <utility>
# include <iostream>
# include <algorithm>
# include <functional>

///
///这是std :: vector的包装。唯一的目的是在
///向量构造在称为。
///它包装了索引访问运算符[]和size()方法,它们分别是 
///对以后的ET实施很重要。
///

// std :: vector包装器。
template<typename ScalarType> class Vector
{
public:
  explicit Vector() { std::cout << "ctor called.\n"; };
  explicit Vector(int size): _vec(size) { std::cout << "ctor called.\n"; };
  explicit Vector(const std::vector<ScalarType> &vec): _vec(vec)
  { std::cout << "ctor called.\n"; };
  
  Vector(const Vector<ScalarType> & vec): _vec{vec()}
  { std::cout << "copy ctor called.\n"; };
  Vector(Vector<ScalarType> && vec): _vec(std::move(vec()))
  { std::cout << "move ctor called.\n"; };

  Vector<ScalarType> & operator=(const Vector<ScalarType> &) = default;
  Vector<ScalarType> & operator=(Vector<ScalarType> &&) = default;

  decltype(auto) operator[](int indx) { return _vec[indx]; }
  decltype(auto) operator[](int indx) const { return _vec[indx]; }

  decltype(auto) operator()() & { return (_vec); };        
  decltype(auto) operator()() const & { return (_vec); };  
  Vector<ScalarType> && operator()() && { return std::move(*this); }

  int size() const { return _vec.size(); }
  
private:
  std::vector<ScalarType> _vec;
};

///
///这些是运算符+(向量加运算)的常规重载
///和运算符*(向量内积运算)而不使用表达式
///模板。它们后来用于基准测试。
///

// +(向量加)运算符。
template<typename ScalarType>
auto operator+(const Vector<ScalarType> &lhs, const Vector<ScalarType> &rhs)
{
  assert(lhs().size() == rhs().size() &&
         "error: ops plus -> lhs and rhs size mismatch.");
  
  std::vector<ScalarType> _vec;
  _vec.resize(lhs().size());
  std::transform(std::cbegin(lhs()), std::cend(lhs()),
                 std::cbegin(rhs()), std::begin(_vec),
                 std::plus<>());
  return Vector<ScalarType>(std::move(_vec));
}

// *(向量内积)运算符。
template<typename ScalarType>
auto operator*(const Vector<ScalarType> &lhs, const Vector<ScalarType> &rhs)
{
  assert(lhs().size() == rhs().size() &&
         "error: ops multiplies -> lhs and rhs size mismatch.");
  
  std::vector<ScalarType> _vec;
  _vec.resize(lhs().size());
  std::transform(std::cbegin(lhs()), std::cend(lhs()),
                 std::cbegin(rhs()), std::begin(_vec),
                 std::multiplies<>());
  return Vector<ScalarType>(std::move(_vec));
}

#endif //!EXPR_VEC


文件expr.hh:用于元素操作(向量加和向量内积)的表达式模板的实现


让我们将其分为几部分。

  1. 第1节为所有表达式实现了一个基类。它采用了好奇重复模板模式(CRTP)。

  2. 第2节实现了第一个PAE:一个terminal,它只是输入数据结构的包装器(const引用),其中包含用于计算的实际输入值。

  3. 第3节实现了第二个PAEbinary_operation,它是稍后用于vector_plus和vector_innerprod的类模板。它是由参数化操作的类型左手端PAE右手边PAE。实际的计算在索引访问运算符中进行编码。

  4. 第4节将vector_plus和vector_innerprod操作定义为逐元素操作。它还使PAE的运算符+和*重载:这样,这两个操作也将返回PAE

#ifndef EXPR_EXPR
# define EXPR_EXPR
      

/// Fwd声明。
template<typename> class Vector;

namespace expr
{


/// -----------------------------------------
///
///第1节。
///
///第一部分是各种表达式的基类模板。它         
/// /使用好奇重复发生的模板模式,该模式可实例化 
///继承自它的任何一种表达式结构。
///
/// -----------------------------------------


  ///所有表达式的基类。
  template<typename Expr> class expr_base
  {
  public:
    const Expr& self() const { return static_cast<const Expr&>(*this); }
    Expr& self() { return static_cast<Expr&>(*this); }

  protected:
    explicit expr_base() {};
    int size() const { return self().size_impl(); }
    auto operator[](int indx) const { return self().at_impl(indx); }
    auto operator()() const { return self()(); };
  };
  

/// -----------------------------------------
///
/// The following section 2 & 3 are abstractions of pure algebraic expressions (PAE).
///可以使用operator()将任何PAE转换为真实对象实例: 
///此转换过程完成了真正的计算。

///
///第2节。终端
///
///终端是将const引用包装到Vector数据的抽象 
/// 结构体。它继承自expr_base,因此提供了一个统一的接口
///将向量包装到PAE中。
///
///它提供size()方法,通过at_impl()进行索引访问和转换
///通过()运算符引用的对象。
/// 
///可能不需要用户定义的数据结构具有终端 
///包装器,因为用户定义的结构可以继承expr_base,因此消除了
///需要提供此类终端包装。 
///
/// -----------------------------------------


  ///底层数据结构的通用包装器。
  template<typename DataType> class terminal: expr_base<terminal<DataType>>
  {
  public:
    using base_type = expr_base<terminal<DataType>>;
    using base_type::size;
    using base_type::operator[];
    friend base_type;
    
    explicit terminal(const DataType &val): _val(val) {}
    int size_impl() const { return _val.size(); };
    auto at_impl(int indx) const { return _val[indx]; };
    decltype(auto) operator()() const { return (_val); }
    
  private:
    const DataType &_val;
  };


/// -----------------------------------------
///
///第3节。二进制运算表达式。
///
///这是任何二进制表达式的PAE抽象。同样,它继承自 
/// expr_base。
///
///它提供size()方法,通过at_impl()进行索引访问和转换
///通过()运算符引用的对象。 Each call to the at_impl() method is
///元素明智的计算。
/// 
/// -----------------------------------------


  ///二进制操作的通用包装器(按元素)。
  template<typename Ops, typename lExpr, typename rExpr>
  class binary_ops: public expr_base<binary_ops<Ops,lExpr,rExpr>>
  {
  public:
    using base_type = expr_base<binary_ops<Ops,lExpr,rExpr>>;
    using base_type::size;
    using base_type::operator[];
    friend base_type;
    
    explicit binary_ops(const Ops &ops, const lExpr &lxpr, const rExpr &rxpr)
      : _ops(ops), _lxpr(lxpr), _rxpr(rxpr) {};
    int size_impl() const { return _lxpr.size(); };

    ///这对索引indx进行逐元素计算。
    auto at_impl(int indx) const { return _ops(_lxpr[indx], _rxpr[indx]); };

    ///从任意expr转换为具体数据类型。它评估
    ///所有索引的逐元素计算。
    template<typename DataType> operator DataType()
    {
      DataType _vec(size());
      for(int _ind = 0; _ind < _vec.size(); ++_ind)
        _vec[_ind] = (*this)[_ind];
      return _vec;
    }
    
  private: /// ops和expr被认为便宜复制。
    Ops   _ops;
    lExpr _lxpr;
    rExpr _rxpr;
  };


/// -----------------------------------------
///第4节。
///
///以下两个结构定义了PAE上的代数运算:这里仅是向量 
/// plus和矢量内积被实现。 
///
///首先,定义一些元素操作:换句话说,vec_plus和 
/// vec_prod作用于Vector中的元素,但不作用于整个Vector。 
///
/// Then, operator + & * are overloaded on PAEs, such that: + & * operations on PAEs         
///还返回PAE。
///
/// -----------------------------------------


  ///按元素加运算。
  struct vec_plus_t
  {
    constexpr explicit vec_plus_t() = default; 
    template<typename LType, typename RType>
    auto operator()(const LType &lhs, const RType &rhs) const
    { return lhs+rhs; }
  };
  
  ///按元素进行内部产品操作。
  struct vec_prod_t
  {
    constexpr explicit vec_prod_t() = default; 
    template<typename LType, typename RType>
    auto operator()(const LType &lhs, const RType &rhs) const
    { return lhs*rhs; }
  };
  
  ///常数加和内部乘积运算符对象。
  constexpr vec_plus_t vec_plus{};
  constexpr vec_prod_t vec_prod{};
  
  ///加运算符对表达式的重载:返回二进制表达式。
  template<typename lExpr, typename rExpr>
  auto operator+(const lExpr &lhs, const rExpr &rhs)
  { return binary_ops<vec_plus_t,lExpr,rExpr>(vec_plus,lhs,rhs); }
  
  ///内部prod运算符对表达式的重载:返回二进制表达式。
  template<typename lExpr, typename rExpr>
  auto operator*(const lExpr &lhs, const rExpr &rhs)
  { return binary_ops<vec_prod_t,lExpr,rExpr>(vec_prod,lhs,rhs); }
  
} //!expr


#endif //!EXPR_EXPR


文件main.cc:测试src文件


# include <chrono>
# include <iomanip>
# include <iostream>
# include "vec.hh"
# include "expr.hh"
# include "boost/core/demangle.hpp"


int main()
{
  using dtype = float;
  constexpr int size = 5e7;
  
  std::vector<dtype> _vec1(size);
  std::vector<dtype> _vec2(size);
  std::vector<dtype> _vec3(size);

  // ...初始化向量的内容。

  Vector<dtype> vec1(std::move(_vec1));
  Vector<dtype> vec2(std::move(_vec2));
  Vector<dtype> vec3(std::move(_vec3));

  unsigned long start_ms_no_ets =
    std::chrono::duration_cast<std::chrono::milliseconds>
    (std::chrono::system_clock::now().time_since_epoch()).count();
  std::cout << "\nNo-ETs evaluation starts.\n";
  
  Vector<dtype> result_no_ets = vec1 + (vec2*vec3);
  
  unsigned long stop_ms_no_ets =
    std::chrono::duration_cast<std::chrono::milliseconds>
    (std::chrono::system_clock::now().time_since_epoch()).count();
  std::cout << std::setprecision(6) << std::fixed
            << "No-ETs. Time eclapses: " << (stop_ms_no_ets-start_ms_no_ets)/1000.0
            << " s.\n" << std::endl;
  
  unsigned long start_ms_ets =
    std::chrono::duration_cast<std::chrono::milliseconds>
    (std::chrono::system_clock::now().time_since_epoch()).count();
  std::cout << "Evaluation using ETs starts.\n";
  
  expr::terminal<Vector<dtype>> vec4(vec1);
  expr::terminal<Vector<dtype>> vec5(vec2);
  expr::terminal<Vector<dtype>> vec6(vec3);
  
  Vector<dtype> result_ets = (vec4 + vec5*vec6);
  
  unsigned long stop_ms_ets =
    std::chrono::duration_cast<std::chrono::milliseconds>
    (std::chrono::system_clock::now().time_since_epoch()).count();
  std::cout << std::setprecision(6) << std::fixed
            << "With ETs. Time eclapses: " << (stop_ms_ets-start_ms_ets)/1000.0
            << " s.\n" << std::endl;
  
  auto ets_ret_type = (vec4 + vec5*vec6);
  std::cout << "\nETs result's type:\n";
  std::cout << boost::core::demangle( typeid(decltype(ets_ret_type)).name() ) << '\n'; 

  return 0;
}

这是-O3 -std=c++14使用GCC 5.3编译时可能的输出:

ctor called.
ctor called.
ctor called.

No-ETs evaluation starts.
ctor called.
ctor called.
No-ETs. Time eclapses: 0.571000 s.

Evaluation using ETs starts.
ctor called.
With ETs. Time eclapses: 0.164000 s.


ETs result's type:
expr::binary_ops<expr::vec_plus_t, expr::terminal<Vector<float> >, expr::binary_ops<expr::vec_prod_t, expr::terminal<Vector<float> >, expr::terminal<Vector<float> > > >

观察结果是:

  • 在这种情况下,使用ET可以显着提高性能(> 3倍)。

  • 消除了创建临时Vector对象的操作。与ET一样,ctor仅被调用一次。

  • Boost :: demangle用于在转换之前可视化ET返回的类型:它清楚地构造了与上面演示的完全相同的表达图。


缺点和警告


  • ET的一个明显缺点是学习曲线,实现的复杂性和代码维护的难度。在上面的示例中,仅考虑了按元素的运算,实现中已经包含了大量的样板,更不用说在现实世界中,在每次计算中都会出现更复杂的代数表达式,并且不再保持按元素的独立性(例如矩阵乘法) ),难度将是指数级的。

  • 使用ET的另一个警告是它们可以很好地与auto关键字一起使用。如上所述,PAE本质上是代理:代理基本上不能与一起使用auto。考虑以下示例:

     auto result = ...;                // 一些昂贵的表达: 
                                      // 自动返回expr图,
                                      // 不是计算值。
    for(auto i = 0; i < 100; ++i)
        ScalrType value = result* ... // 其他一些使用结果的昂贵计算。

在此,在for循环的每次迭代中,将重新评估结果,因为表达式图而不是计算值将传递给for循环。


现有实现ET的


  • boost :: proto是一个功能强大的库,可让您为自己的表达式定义自己的规则和语法,并使用ET执行。

  • Eigen是用于线性代数的库,可使用ET有效地实现各种代数计算。