Boost.Proto: возможно ли для преобразований Proto оценить смешанное выражение матрицы и вектора?

Теперь я пытаюсь научить линейную алгебру компилятора g++, чтобы g++ мог переписать такое выражение, как (matrix * vector)(index) как цикл для оценки выражения. В основном это то, что я ожидаю в качестве следующей статьи последней статьи в серии " Expressive C++". В этой последней статье объясняется, как создать EDSL для добавления векторов, так что я написал еще один EDSL для умножения матрицы на вектор.

Но макрос BOOST_PROTO_DEFINE_OPERATORS не может быть скомпилирован, когда имя домена Proto для моей собственной матрицы и векторных классов передается в качестве первого параметра макроса.

Поэтому мне интересно, возможно ли для преобразований Прото оценить смешанное выражение матричных и векторных объектов. Кажется, нет такого примера кода, который можно скомпилировать, и пример кода в "Адаптации существующих типов к Proto" в руководстве пользователя Proto 1.57.0 неполон, хотя он посвящен тому, как адаптировать существующие матрицы и векторные типы к Proto.

Я в растерянности.. Не могли бы вы дать мне несколько советов или подсказок?

А вот и мой исходный код. Я был бы очень признателен, если бы вы посоветовали мне, как это исправить:

#include <iostream>
#include <boost/proto/proto.hpp>

namespace mpl = boost::mpl;
namespace proto = boost::proto;

namespace LinAlg {
    class Vector;
    class Matrix;

    // Functor for lazy evaluation
    struct ElmOfMatVecMult;

    // The grammar for an expression like ( matrix * vector )(index)
    struct MatVecMultElmGrammar : proto::or_<
        proto::when< proto::multiplies< Matrix, Vector>,
                    proto::_make_function( ElmOfMatVecMult,
                                            proto::_left, proto::_right,
                                            proto::_state) >
    > {};

    // The grammar for a vector expression
    // ( Now I consider just one form : matrix * vector . )
    struct VecExprGrammar : proto::or_<
        proto::when< proto::function< MatVecMultElmGrammar, proto::_>,
                    MatVecMultElmGrammar( proto::_left, proto::_right) >,
        proto::multiplies< Matrix, Vector>
    > {};

    template<typename E> struct Expr;

    // The above grammar is associated with this domain.
    struct Domain
        : proto::domain<proto::generator<Expr>, VecExprGrammar>
    {};

    // A wrapper template for linear algebraic expressions
    // including matrices and vectors
    template<typename ExprType>
    struct Expr
        : proto::extends<ExprType, Expr<ExprType>, Domain>
    {
        explicit Expr(const ExprType& e)
            : proto::extends<ExprType, Expr<ExprType>, Domain>(e)
        {}
    };

    // Testing if data in an heap array can be a vector object
    class Vector {
        private:
            int sz;
            double* data;

    public:
        template <typename Sig> struct result;

        template <typename This, typename T>
        struct result< This(T) > { typedef double type; };

        typedef double result_type;

        explicit Vector(int sz_ = 1, double iniVal = 0.0) :
            sz( sz_), data( new double[sz] ) {
            for (int i = 0; i < sz; i++) data[i] = iniVal;
            std::cout << "Created" << std::endl;
        }
        Vector(const Vector& vec) :
            sz( vec.sz), data( new double[sz] ) {
            for (int i = 0; i < sz; i++) data[i] = vec.data[i];
            std::cout << "Copied! " << std::endl;
        }

        ~Vector() {
            delete [] data;
            std::cout << "Deleted" << std::endl;
        }

        // accessing to an element of this vector
        double& operator()(int i) { return data[i]; }
        const double& operator()(int i) const { return data[i]; }
    };


    // Testing if data in an heap array can be a matrix object
    class Matrix
    {
    private:
        int rowSz, colSz;
        double* data;
        double** m;

    public:
        template <typename Signature> struct result;

        template <typename This, typename T>
        struct result< This(T,T) > { typedef double type; };

        explicit Matrix(int rowSize = 1, int columnSize =1,
                    double iniVal = 0.0) :
            rowSz( rowSize), colSz(columnSize),
            data( new double[rowSz*colSz] ), m( new double*[rowSz])
        {
            for (int i = 0; i < rowSz; i++) m[i] = data + i*colSz;
            for (int ri = 0; ri < rowSz; ri++)
                for (int ci = 0; ci < colSz; ci++) m[ri][ci] = iniVal;
            std::cout << "Created" << std::endl;
        }

        Matrix(const Matrix& mat) :
            rowSz( mat.rowSz), colSz( mat.colSz),
            data( new double[rowSz*colSz] ), m( new double*[rowSz])
        {
            for (int i = 0; i < rowSz; i++) m[i] = data + i*colSz;
                for (int ri = 0; ri < rowSz; ri++)
                    for (int ci = 0; ci < colSz; ci++)
                        m[ri][ci] = mat.m[ri][ci];
                std::cout << "Copied! " << std::endl;
        }

        ~Matrix()
        {
            delete [] m;
            delete [] data;
            std::cout << "Deleted" << std::endl;
        }

        int rowSize() const { return rowSz; }
        int columnSize() const { return colSz; }

        // accesing to a vector element
        double& operator()(int ri, int ci) { return m[ri][ci]; }
        const double& operator()(int ri, int ci) const { return m[ri][ci]; }
    };

    // An expression like ( matrix * vector )(index) is transformed
    // into the loop for calculating the dot product between
    // the vector and matrix.
    struct ElmOfMatVecMult
    {
        double operator()( Matrix const& mat, Vector const& vec,
                        int index) const
        {
            double elm = 0.0;
            for (int ci =0;  ci < mat.columnSize(); ci++)
                elm += mat(index, ci) * vec(ci);
            return elm;
        }
    };

    // Define a trait for detecting linear algebraic terminals, to be used
    // by the BOOST_PROTO_DEFINE_OPERATORS macro below.
    template<typename> struct IsExpr  : mpl::false_ {};
    template<> struct IsExpr< Vector> : mpl::true_  {};
    template<> struct IsExpr< Matrix> : mpl::true_  {};

    // This defines all the overloads to make expressions involving
    // Vector and Matrix objects to build Proto's expression templates.
    BOOST_PROTO_DEFINE_OPERATORS(IsExpr, Domain)
}

int main()
{
    using namespace LinAlg;

    proto::_default<> trans;

    Matrix mat( 3, 3);
    Vector vec1(3);

    mat(0,0) = 1.00; mat(0,1) = 1.01; mat(0,2) = 1.02;
    mat(1,0) = 1.10; mat(1,1) = 1.11; mat(1,2) = 1.12;
    mat(2,0) = 1.20; mat(2,1) = 1.21; mat(2,2) = 1.22;

    vec1(0) = 1.0;
    vec1(1) = 2.0;
    vec1(2) = 3.0;

    proto::display_expr( ( mat * vec1)(2) );
    proto::display_expr( VecExprGrammar()( ( mat * vec1)(2) );
    double vecElm = trans( VecExprGrammar()( ( mat * vec1)(2) );

    return 0;
}

1 ответ

Решение

Я как-то подтвердил, что ответ на мой вопрос "да". Хитрость моего исходного кода в том, чтобы определить активную грамматику для замены вектора матрицы * на MatVecMult объект и определить, что MatVecMult как proto::callable объект для создания ленивого оценочного функтора, LazyMatVecMult, Когда это LazyMatVecMult Экземпляр вызывается с operator()(int index), он вычисляет скалярное произведение индексной строки матрицы и вектора.

(Извините за то, что сам ответил на мой вопрос. Но я думаю, что я не одинок, чтобы наткнуться на этот камень о том, как сделать основанный на Proto EDSL для векторных и матричных выражений.)

Этот исходный код подтвержден для работы с g++ 4.8.4:

#include <iostream>
#include <boost/proto/proto.hpp>

namespace mpl = boost::mpl;
namespace proto = boost::proto;

namespace LinAlg {
    class Vector;
    class Matrix;


    // Callable transform object to make a proto exression
    // for lazily evaluationg a. multiplication
    struct MatVecMult;

    // The grammar for the multiplication of a matrix and a vector
    struct MatVecMultGrammar : proto::or_<
        proto::when<
            proto::multiplies< proto::terminal< Matrix> ,
                                proto::terminal< Vector> >,
            MatVecMult( proto::_value( proto::_left),
                        proto::_value( proto::_right) )
        >
    > {};

    // The grammar for a vector expression
    // ( Now I consider just one form : matrix * vector . )
    struct VecExprGrammar : proto::or_<
        proto::function< MatVecMultGrammar, proto::_>,
        MatVecMultGrammar
    > {};


    // A wrapper for a linear algebraic expression
    template<typename E> struct Expr;

    // The above grammar is associated with this domain.
    struct Domain
        : proto::domain<proto::generator<Expr>, VecExprGrammar>
    {};

    // A wrapper template for linear algebraic expressions
    // including matrices and vectors
    template<typename ExprType>
    struct Expr
        : proto::extends<ExprType, Expr<ExprType>, Domain>
    {
        /* typedef double result_type; */

        explicit Expr(const ExprType& e)
            : proto::extends<ExprType, Expr<ExprType>, Domain>(e)
        {}
    };


    // Testing if data in an heap array can be a vector object
    class Vector {
        private:
            int sz;
            double* data;

    public:
        explicit Vector(int sz_ = 1, double iniVal = 0.0) :
            sz( sz_), data( new double[sz] ) {
            for (int i = 0; i < sz; i++) data[i] = iniVal;
            std::cout << "Created" << std::endl;
        }
        Vector(const Vector& vec) :
            sz( vec.sz), data( new double[sz] ) {
            for (int i = 0; i < sz; i++) data[i] = vec.data[i];
            std::cout << "Copied! " << std::endl;
        }

        ~Vector() {
            delete [] data;
            std::cout << "Deleted" << std::endl;
        }

        // accessing to an element of this vector
        double& operator()(int i) { return data[i]; }
        const double& operator()(int i) const { return data[i]; }
    };


    // Testing if data in an heap array can be a matrix object
    class Matrix
    {
    private:
        int rowSz, colSz;
        double* data;
        double** m;

    public:

        explicit Matrix(int rowSize = 1, int columnSize =1,
                    double iniVal = 0.0) :
            rowSz( rowSize), colSz(columnSize),
            data( new double[rowSz*colSz] ), m( new double*[rowSz])
        {
            for (int i = 0; i < rowSz; i++) m[i] = data + i*colSz;
            for (int ri = 0; ri < rowSz; ri++)
                for (int ci = 0; ci < colSz; ci++) m[ri][ci] = iniVal;
            std::cout << "Created" << std::endl;
        }

        Matrix(const Matrix& mat) :
            rowSz( mat.rowSz), colSz( mat.colSz),
            data( new double[rowSz*colSz] ), m( new double*[rowSz])
        {
            for (int i = 0; i < rowSz; i++) m[i] = data + i*colSz;
                for (int ri = 0; ri < rowSz; ri++)
                    for (int ci = 0; ci < colSz; ci++)
                        m[ri][ci] = mat.m[ri][ci];
                std::cout << "Copied! " << std::endl;
        }

        ~Matrix()
        {
            delete [] m;
            delete [] data;
            std::cout << "Deleted" << std::endl;
        }

        int rowSize() const { return rowSz; }
        int columnSize() const { return colSz; }

        // accesing to a vector element
        double& operator()(int ri, int ci) { return m[ri][ci]; }
        const double& operator()(int ri, int ci) const { return m[ri][ci]; }

    };

    // Lazy function object for evaluating an element of
    // the resultant vector from the multiplication of
    // a matrix and vector objects.
    //
    // An expression like ( matrix * vector )(index) is transformed
    // into the loop for calculating the dot product between
    // the index'th row of the matrix and the vector.
    struct LazyMatVecMult
    {
        Matrix const& m;
        Vector const& v;
        int mColSz;

        typedef double result_type;
        // typedef mpl::int_<1> proto_arity;

        explicit LazyMatVecMult(Matrix const& mat, Vector const& vec) :
        m( mat), v( vec), mColSz(mat.rowSize()) {}

        LazyMatVecMult(LazyMatVecMult const& lazy) :
            m(lazy.m), v(lazy.v), mColSz(lazy.mColSz) {}

        result_type operator()(int index) const
        {
            result_type elm = 0.0;
            for (int ci =0;  ci < mColSz; ci++)
                elm += m(index, ci) * v(ci);
            return elm;
        }
    };

    // Callable transform object to make the lazy functor
    // a proto exression for lazily evaluationg the multiplication
    // of a matrix and a vector .
    struct MatVecMult : proto::callable
    {
        typedef proto::terminal< LazyMatVecMult >::type result_type;

        result_type
        operator()( Matrix const& mat, Vector const& vec) const
        {
            return proto::as_expr( LazyMatVecMult(mat, vec) );
        }
    };

    // Define a trait for detecting linear algebraic terminals, to be used
    // by the BOOST_PROTO_DEFINE_OPERATORS macro below.
    template<typename> struct IsExpr  : mpl::false_ {};
    template<> struct IsExpr< Vector> : mpl::true_  {};
    template<> struct IsExpr< Matrix> : mpl::true_  {};
    // template<> struct IsExpr< MatVecMult> : mpl::true_  {};
    template<> struct IsExpr< LazyMatVecMult> : mpl::true_  {};

    // This defines all the overloads to make expressions involving
    // Vector and Matrix objects to build Proto's expression templates.
    BOOST_PROTO_DEFINE_OPERATORS(IsExpr, Domain)
    // BOOST_PROTO_DEFINE_OPERATORS(IsExpr, proto::default_domain)


}



int main()
{
    using namespace LinAlg;

    proto::_default<> trans;

    Matrix mat( 3, 3);
    Vector vec1(3), vec2(3);

    mat(0,0) = 1.00; mat(0,1) = 1.01; mat(0,2) = 1.02;
    mat(1,0) = 1.10; mat(1,1) = 1.11; mat(1,2) = 1.12;
    mat(2,0) = 1.20; mat(2,1) = 1.21; mat(2,2) = 1.22;

    vec1(0) = 1.0;
    vec1(1) = 2.0;
    vec1(2) = 3.0;

    std::cout << " mat * vec1" << std::endl;
    proto::display_expr( mat * vec1 );
    std::cout << " MatVecMultGrammar()( mat * vec1) " << std::endl;
    proto::display_expr( MatVecMultGrammar()( mat * vec1) );

    std::cout << "( mat * vec1)(2)" << std::endl;
    proto::display_expr( ( mat * vec1)(2) );
    std::cout << "VecExprGrammar()( ( mat * vec1)(2) )" << std::endl;
    proto::display_expr( VecExprGrammar()( ( mat * vec1)(2) ) );
    double elm2 = trans( VecExprGrammar()( ( mat * vec1)(2) ) );

    std::cout << "( mat * vec1)(2) = " << elm2 << std::endl;

    return 0;
}
Другие вопросы по тегам