/****************************************************************************** * * AMDiS - Adaptive multidimensional simulations * * Copyright (C) 2013 Dresden University of Technology. All Rights Reserved. * Web: https://fusionforge.zih.tu-dresden.de/projects/amdis * * Authors: * Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al. * * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. * * * This file is part of AMDiS * * See also license.opensource.txt in the distribution. * ******************************************************************************/ /** \file vec_functors.hpp */ #ifndef AMDIS_VEC_FUNCTORS_HPP #define AMDIS_VEC_FUNCTORS_HPP #include "functor_expr.hpp" #include "operations/norm.hpp" #include "operations/product.hpp" /** * This file provides expressions for vectors and matrices: * (where v_expr is an expression of vector type, and m_expr an * expression of matrix type) * * two_norm(v_expr) ... the 2-norm of a vector v: result = sqrt(v^H * v) * one_norm(v_expr) ... the 1-norm of a vector v: result = sum_i(abs(v_i)) * one_norm(m_expr) ... the 1-norm of a matrix m: result = max_j(sum_i(abs(m_ij)) * p_norm

(v_expr) .. the P-norm of a vector v: result = [sum_i(abs(v_i)^P)]^(1/P) * **/ namespace AMDiS { namespace traits { template struct is_always_true : boost::mpl::true_ {}; } namespace result_of { // helper class for UnaryExpr template class Functor, typename Term, template class Condition1 = traits::is_always_true, template class Condition2 = traits::is_always_true, template class Condition3 = traits::is_always_true> struct UnaryExpr : boost::enable_if < typename boost::mpl::and_ < typename traits::is_expr::type, typename boost::mpl::or_ < typename Condition1::type, typename Condition2::type, typename Condition3::type >::type >, expressions::Function1, Term> > {}; // helper class for UnaryExpr template class Condition1 = traits::is_always_true, template class Condition2 = traits::is_always_true, template class Condition3 = traits::is_always_true> struct UnaryExprFull : boost::enable_if < typename boost::mpl::and_ < typename traits::is_expr::type, typename boost::mpl::or_ < typename Condition1::type, typename Condition2::type, typename Condition3::type >::type >, expressions::Function1 > {}; // helper class for BinaryExpr template class Functor, typename Term1, typename Term2, template class Condition1 = traits::is_always_true, template class Condition2 = traits::is_always_true, template class Condition3 = traits::is_always_true, typename Expr1 = typename traits::to_expr::type, typename Expr2 = typename traits::to_expr::type> struct BinaryExpr : boost::enable_if < typename boost::mpl::and_ < typename traits::is_valid_arg2::type, typename boost::mpl::or_ < typename Condition1::type, typename Condition2::type, typename Condition3::type >::type, typename boost::mpl::or_ < typename Condition1::type, typename Condition2::type, typename Condition3::type >::type >::type, expressions::Function2 < Functor, Expr1, Expr2 > > {}; } // end namespace result_of /// expression with one vector template class Functor, typename Term> inline typename result_of::UnaryExpr::type unary_expr(const Term& t) { typedef Functor< typename Term::value_type > F; return function_(F(), t); } template inline typename result_of::UnaryExprFull::type unary_expr_full(const Term& t) { return function_(F(), t); } /// expression with two vectors v1, v2 template class Functor, typename Term1, typename Term2> inline typename result_of::BinaryExpr::type binary_expr(const Term1& t1, const Term2& t2) { typedef typename traits::to_expr::to Expr1; typedef typename traits::to_expr::to Expr2; typedef Functor< typename Expr1::type::value_type, typename Expr2::type::value_type > F; return function_(F(), Expr1::get(t1), Expr2::get(t2)); } // two_norm of vectors // _____________________________________________________________________________ namespace result_of { template struct two_norm : result_of::UnaryExpr {}; } /// the 2-norm of a vector v: result = sqrt(v^H * v) template inline typename result_of::two_norm::type two_norm_dispatch(const Term& t, tag::expression) { return unary_expr< functors::TwoNorm >(t); } // one_norm of vectors and matrices // _____________________________________________________________________________ namespace result_of { template struct one_norm : result_of::UnaryExpr {}; } /// the 1-norm of a vector v: result = max_i(abs(v_i)) template inline typename result_of::one_norm::type one_norm_dispatch(const Term& t, tag::expression) { return unary_expr< functors::OneNorm >(t); } // p_norm

of vectors // _____________________________________________________________________________ namespace result_of { template struct p_norm : result_of::UnaryExprFull, T, traits::is_vector> {}; } /// the 2-norm of a vector v: result = sqrt(v^H * v) template inline typename result_of::p_norm::type p_norm_dispatch(const Term& t, tag::expression) { return unary_expr_full< functors::PNorm >(t); } // unary dot (v' * v) of a vector v // _____________________________________________________________________________ namespace result_of { template struct unary_dot : result_of::UnaryExpr {}; } /// the 2-norm squared of a vector v: result = v^H * v template inline typename result_of::unary_dot::type unary_dot_dispatch(const Term& t, tag::expression) { return unary_expr< functors::UnaryDot >(t); } // inner product of two vectors // _____________________________________________________________________________ namespace result_of { template struct dot : result_of::BinaryExpr {}; template struct dot : result_of::BinaryExpr {}; template struct dot : result_of::BinaryExpr {}; } /// inner/dot product of two vectors v1, v2: result = v1^T * v2 template inline typename boost::enable_if< typename traits::is_expr::type, typename result_of::dot::type >::type dot_dispatch(const Term1& t1, const Term2& t2, tag::expression tag1_, Tag2 tag2_) { return binary_expr(t1, t2); } template inline typename boost::enable_if_c< (traits::is_expr::type::value && !(traits::is_expr::type::value)), typename result_of::dot::type >::type dot_dispatch(const Term1& t1, const Term2& t2, Tag1 tag1_, tag::expression tag2_) { return binary_expr(t1, t2); } // cross-product of two vectors // _____________________________________________________________________________ namespace result_of { template struct cross : result_of::BinaryExpr {}; template struct cross : result_of::BinaryExpr {}; template struct cross : result_of::BinaryExpr {}; } /// cross product of two vectors v1, v2: result = v1 x v2 template inline typename boost::enable_if< typename traits::is_expr::type, typename result_of::cross::type >::type cross_dispatch(const Term1& t1, const Term2& t2, tag::expression tag1_, Tag2 tag2_) { return binary_expr(t1, t2); } template inline typename boost::enable_if_c< (traits::is_expr::type::value && !(traits::is_expr::type::value)), typename result_of::cross::type >::type cross_dispatch(const Term1& t1, const Term2& t2, Tag1 tag1_, tag::expression tag2_) { return binary_expr(t1, t2); } // generator for a diagonal matrix from a vector // _____________________________________________________________________________ namespace expressions { template struct DiagonalMat : public FunctorBase {}; template struct DiagonalMat > : public FunctorBase { typedef WorldMatrix result_type; int getDegree(int d0) const { return d0; } result_type operator()(const WorldVector &v) const { T zero; nullify(zero); result_type matrix(DEFAULT_VALUE, zero); for (int i = 0; i < v.getSize(); ++i) matrix[i][i] = v[i]; return matrix; } }; template struct DiagonalMat > : public FunctorBase { typedef DimMat result_type; int getDegree(int d0) const { return d0; } result_type operator()(const DimVec &v) const { T zero; nullify(zero); result_type matrix(v.getDim(), DEFAULT_VALUE, zero); for (int i = 0; i < v.getSize(); ++i) matrix[i][i] = v[i]; return matrix; } }; template struct DiagonalMat > : public FunctorBase { typedef Matrix result_type; int getDegree(int d0) const { return d0; } result_type operator()(const DimVec &v) const { T zero; nullify(zero); result_type matrix(v.getSize(), v.getSize()); matrix.set(zero); for (int i = 0; i < v.getSize(); ++i) matrix[i][i] = v[i]; return matrix; } }; template struct DiagonalMat > : public FunctorBase { typedef mtl::matrix::compressed2D > result_type; int getDegree(int d0) const { return d0; } template result_type operator()(const Vector &v) const { return diagonal(v); } }; } /// create diagonal matrix from vector template inline typename result_of::UnaryExpr::type diagonal(const Term& t) { return unary_expr(t); } // outer product / dyadic product of two vectors to create a matrix // _____________________________________________________________________________ /// outer product of two vectors v1, v2: result = v1 * v2^T template inline typename result_of::BinaryExpr::type outer(const Term1& t1, const Term2& t2) { return binary_expr(t1, t2); } // extract a component of a vector/matrix // _____________________________________________________________________________ namespace expressions { template struct VecComponent : public FunctorBase { typedef typename traits::category::value_type result_type; VecComponent(int I_) : I(I_) {} int getDegree(int d0) const { return d0; } result_type operator()(const T &v) const { return v[I]; } private: int I; }; template struct MatComponent : public FunctorBase { typedef typename traits::category::value_type result_type; MatComponent(int I_, int J_) : I(I_), J(J_) {} int getDegree(int d0) const { return d0; } result_type operator()(const T &m) const { return m[I][J]; } private: int I; int J; }; } template typename result_of::UnaryExpr::type at(const Term& t, int I) { return function_(expressions::VecComponent(I), t); } template typename result_of::UnaryExpr::type at(const Term& t, int I, int J) { return function_(expressions::MatComponent(I, J), t); } // transpose a matrix // _____________________________________________________________________________ namespace expressions { template struct MatTranspose : public FunctorBase { typedef Mat result_type; int getDegree(int d0) const { return d0; } result_type operator()(const Mat &m) const { Mat result; for (size_t r = 0; r < num_rows(m); r++) for (size_t c = 0; c < num_cols(m); c++) result[c][r] = m[r][c]; return result; } }; } template typename result_of::UnaryExpr::type trans(const Term& t) { return function_(expressions::MatTranspose(), t); } } #endif // AMDIS_VEC_FUNCTORS_HPP