From ca89f782e3b39699e2dfb6b5cb87e90b025b8372 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Sat, 14 Dec 2019 20:15:27 +0800 Subject: [PATCH 01/22] Initial implementation --- stan/math/prim/mat/err/check_column_index.hpp | 4 +- stan/math/prim/mat/err/check_row_index.hpp | 4 +- stan/math/prim/mat/fun/head.hpp | 25 ++- stan/math/prim/mat/fun/log_softmax.hpp | 12 +- stan/math/prim/mat/fun/log_sum_exp.hpp | 17 +++ .../prim/mat/vectorize/apply_vector_unary.hpp | 144 ++++++++++++++++++ test/unit/math/prim/mat/fun/head_test.cpp | 7 + .../math/prim/mat/fun/log_softmax_test.cpp | 25 ++- .../math/prim/mat/fun/log_sum_exp_test.cpp | 27 +++- 9 files changed, 246 insertions(+), 19 deletions(-) create mode 100644 stan/math/prim/mat/vectorize/apply_vector_unary.hpp diff --git a/stan/math/prim/mat/err/check_column_index.hpp b/stan/math/prim/mat/err/check_column_index.hpp index 13a8127765d..41d04212b97 100644 --- a/stan/math/prim/mat/err/check_column_index.hpp +++ b/stan/math/prim/mat/err/check_column_index.hpp @@ -26,9 +26,9 @@ namespace math { * @param i column index to check * @throw std::out_of_range if index is an invalid column */ -template +template inline void check_column_index(const char* function, const char* name, - const Eigen::Matrix& y, size_t i) { + const Eigen::MatrixBase& y, size_t i) { if (i >= stan::error_index::value && i < static_cast(y.cols()) + stan::error_index::value) { return; diff --git a/stan/math/prim/mat/err/check_row_index.hpp b/stan/math/prim/mat/err/check_row_index.hpp index 7df742b7cc7..c0d2ec006f1 100644 --- a/stan/math/prim/mat/err/check_row_index.hpp +++ b/stan/math/prim/mat/err/check_row_index.hpp @@ -23,9 +23,9 @@ namespace math { * @param i row index to check * @throw std::out_of_range if the index is out of range. */ -template +template inline void check_row_index(const char* function, const char* name, - const Eigen::Matrix& y, size_t i) { + const Eigen::MatrixBase& y, size_t i) { if (i >= stan::error_index::value && i < static_cast(y.rows()) + stan::error_index::value) { return; diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index 8138793a467..9d077bdd89e 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -5,6 +5,7 @@ #include #include #include +#include #include namespace stan { @@ -19,7 +20,7 @@ namespace math { * @param n Size of return. * @return The first n elements of v. * @throw std::out_of_range if n is out of range. - */ + *//* template inline Eigen::Matrix head( const Eigen::Matrix& v, size_t n) { @@ -27,6 +28,20 @@ inline Eigen::Matrix head( check_row_index("head", "n", v, n); } return v.head(n); +}*/ + +template +inline auto head(const T& x, const T2& n) { + return apply_vector_unary::apply_scalar(x, n, [](auto& v, auto& m){ + if (m != 0){ + if (v.rows() == 1){ + check_column_index("head", "n", v, m); + } else { + check_row_index("head", "n", v, m); + } + } + return v.head(m); + }); } /** @@ -38,7 +53,7 @@ inline Eigen::Matrix head( * @param n Size of return row vector. * @return The first n elements of rv. * @throw std::out_of_range if n is out of range. - */ + *//* template inline Eigen::Matrix head( const Eigen::Matrix& rv, size_t n) { @@ -46,7 +61,7 @@ inline Eigen::Matrix head( check_column_index("head", "n", rv, n); } return rv.head(n); -} +}*/ /** * Return the specified number of elements as a standard vector @@ -57,7 +72,7 @@ inline Eigen::Matrix head( * @param n Size of return. * @return The first n elements of sv. * @throw std::out_of_range if n is out of range. - */ + *//* template std::vector head(const std::vector& sv, size_t n) { if (n != 0) { @@ -69,7 +84,7 @@ std::vector head(const std::vector& sv, size_t n) { s.push_back(sv[i]); } return s; -} +}*/ } // namespace math } // namespace stan diff --git a/stan/math/prim/mat/fun/log_softmax.hpp b/stan/math/prim/mat/fun/log_softmax.hpp index 010c78773e3..9dfc390ed71 100644 --- a/stan/math/prim/mat/fun/log_softmax.hpp +++ b/stan/math/prim/mat/fun/log_softmax.hpp @@ -4,6 +4,8 @@ #include #include #include +#include + namespace stan { namespace math { @@ -37,12 +39,12 @@ namespace math { * @return Unit simplex result of the softmax transform of the vector. */ template -inline Eigen::Matrix log_softmax( - const Eigen::Matrix& v) { - check_nonzero_size("log_softmax", "v", v); - return v.array() - log_sum_exp(v); +inline auto log_softmax(const T& x) { + return apply_vector_unary::apply(x, [](auto& v){ + check_nonzero_size("log_softmax", "v", v); + return v.array() - log_sum_exp(v); + }); } - } // namespace math } // namespace stan #endif diff --git a/stan/math/prim/mat/fun/log_sum_exp.hpp b/stan/math/prim/mat/fun/log_sum_exp.hpp index e0c28f32df1..07ed5ad1910 100644 --- a/stan/math/prim/mat/fun/log_sum_exp.hpp +++ b/stan/math/prim/mat/fun/log_sum_exp.hpp @@ -2,6 +2,7 @@ #define STAN_MATH_PRIM_MAT_FUN_LOG_SUM_EXP_HPP #include +#include #include #include #include @@ -23,6 +24,7 @@ namespace math { * @param[in] x Matrix of specified values * @return The log of the sum of the exponentiated vector values. */ + /* template double log_sum_exp(const Eigen::Matrix& x) { if (x.size() == 0) { @@ -34,6 +36,21 @@ double log_sum_exp(const Eigen::Matrix& x) { return max; } return max + std::log((x.array() - max).exp().sum()); +}*/ + +template +inline auto log_sum_exp(const T& x) { + return apply_vector_unary::reduce(x, [](auto& v){ + if (v.size() == 0) { + return -std::numeric_limits::infinity(); + } + + const double max = v.maxCoeff(); + if (!std::isfinite(max)) { + return max; + } + return max + std::log((v.array() - max).exp().sum()); + }); } } // namespace math diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp new file mode 100644 index 00000000000..df6c5e17c1b --- /dev/null +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -0,0 +1,144 @@ +#ifndef STAN_MATH_PRIM_MAT_VECTORIZE_APPLY_VECTOR_UNARY_HPP +#define STAN_MATH_PRIM_MAT_VECTORIZE_APPLY_VECTOR_UNARY_HPP + +#include +#include +#include + +namespace stan { +namespace math { + +template +struct apply_vector_unary {}; + +template +struct apply_vector_unary> { + template + static inline auto apply(const T& x, const F& f) { + return f(x); + } + + template + static inline auto apply_scalar(const T& x, const T2& y, const F& f) { + return f(x, y); + } + + template + static inline auto reduce(const T& x, const F& f) { + return f(x); + } +}; + +template +struct apply_vector_unary> { + using scalar_t = scalar_type_t; + using map_t = + typename Eigen::Map>; + + template + static inline std::vector apply(const T& x, const F& f) { + Eigen::Matrix result + = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); + return std::vector(result.data(), + result.data() + result.size()); + } + + template + static inline std::vector apply_scalar(const T& x, const T2& y, const F& f) { + Eigen::Matrix result + = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), y, f); + return std::vector(result.data(), + result.data() + result.size()); + } + + template + static inline scalar_t reduce(const T& x, const F& f) { + return apply_vector_unary::reduce(as_column_vector_or_scalar(x), f); + } +}; + +template +struct apply_vector_unary> { + using eigen_t = typename T::value_type; + using scalar_t = typename eigen_t::Scalar; + using return_t = std::vector>; + + template + static inline return_t apply(const T& x, const F& f) { + size_t x_size = x.size(); + return_t result(x_size); + for(size_t i = 0; i < x_size; ++i) + result[i] = apply_vector_unary::apply(x[i], f); + return result; + } + + template + static inline return_t apply_scalar(const T& x, const T2& y, const F& f) { + scalar_seq_view y_vec(y); + size_t x_size = x.size(); + return_t result(x_size); + for(size_t i = 0; i < x_size; ++i) + result[i] = apply_vector_unary::apply_scalar(x[i], y_vec[i], f); + return result; + } + + template + static inline std::vector reduce(const T& x, const F& f) { + size_t x_size = x.size(); + std::vector result(x_size); + for(size_t i = 0; i < x_size; ++i) + result[i] = apply_vector_unary::reduce(x[i], f); + return result; + } +}; + +template +struct apply_vector_unary> { + using scalar_t = scalar_type_t; + using return_t = typename std::vector>; + using map_t = + typename Eigen::Map>; + + template + static inline return_t apply(const T& x, const F& f) { + size_t x_size = x.size(); + return_t result(x_size); + Eigen::Matrix inter; + for(size_t i = 0; i < x_size; ++i){ + inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), f); + result[i] = std::vector(inter.data(), + inter.data() + inter.size()); + } + return result; + } + + template + static inline return_t apply_scalar(const T& x, const T2& y, const F& f) { + scalar_seq_view y_vec(y); + size_t x_size = x.size(); + return_t result(x_size); + Eigen::Matrix inter; + for(size_t i = 0; i < x_size; ++i){ + inter = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x[i]), y_vec[i], f); + result[i] = std::vector(inter.data(), + inter.data() + inter.size()); + } + return result; + } + + template + static inline std::vector reduce(const T& x, const F& f) { + size_t x_size = x.size(); + std::vector result(x_size); + for(size_t i = 0; i < x_size; ++i) + result[i] = apply_vector_unary::reduce(as_column_vector_or_scalar(x[i]), f); + + return result; + } +}; + +} // namespace math +} // namespace stan +#endif diff --git a/test/unit/math/prim/mat/fun/head_test.cpp b/test/unit/math/prim/mat/fun/head_test.cpp index 6f5910e3291..53b56baa817 100644 --- a/test/unit/math/prim/mat/fun/head_test.cpp +++ b/test/unit/math/prim/mat/fun/head_test.cpp @@ -8,6 +8,7 @@ TEST(MathMatrixHead, HeadVector1) { v << 1, 2, 3; EXPECT_EQ(0, head(v, 0).size()); } + TEST(MathMatrixHead, HeadVector2) { using stan::math::head; Eigen::VectorXd v(3); @@ -53,6 +54,10 @@ TEST(MathMatrixHead, HeadRowVector4) { using stan::math::head; Eigen::RowVectorXd v(3); v << 1, 2, 3; + std::vector vind{1,2,1}; + + std::vector st_v{v, v, v}; + std::vector st_t = head(st_v, vind); Eigen::RowVectorXd v01 = head(v, 2); EXPECT_EQ(2, v01.size()); @@ -82,6 +87,8 @@ TEST(MathMatrixHead, HeadStdVector3) { v.push_back(1); v.push_back(2); v.push_back(3); + std::vector> st_v{v, v, v}; + std::vector> st_t = head(st_v, 2); EXPECT_THROW(head(v, 4), std::out_of_range); } TEST(MathMatrixHead, HeadStdVector4) { diff --git a/test/unit/math/prim/mat/fun/log_softmax_test.cpp b/test/unit/math/prim/mat/fun/log_softmax_test.cpp index 6a5ef9dafd8..2d28873958b 100644 --- a/test/unit/math/prim/mat/fun/log_softmax_test.cpp +++ b/test/unit/math/prim/mat/fun/log_softmax_test.cpp @@ -33,9 +33,32 @@ TEST(MathMatrixPrimMat, log_softmax) { // x << 0.0; // test_log_softmax(x); + std::vector in{-1,1}; + std::vector out = log_softmax(in); + stan::math::vector_d x2(2); x2 << -1.0, 1.0; - test_log_softmax(x2); + stan::math::matrix_d m2(2,2); + m2 << -1.0, 1.0, -1.0, 1.0; + stan::math::vector_d x2_out = log_softmax(x2); + + EXPECT_FLOAT_EQ(out[0],x2_out[0]); + EXPECT_FLOAT_EQ(out[1],x2_out[1]); + + x2_out = log_softmax(m2.diagonal()); + + EXPECT_FLOAT_EQ(out[0],x2_out[0]); + EXPECT_FLOAT_EQ(out[1],x2_out[1]); + + std::vector invec{x2, x2}; + std::vector outvec = log_softmax(invec); + std::vector> instvec{in, in}; + std::vector> outstvec = log_softmax(instvec); + + EXPECT_FLOAT_EQ(outvec[0][0],outstvec[0][0]); + EXPECT_FLOAT_EQ(outvec[0][1],outstvec[0][1]); + EXPECT_FLOAT_EQ(outvec[1][0],outstvec[1][0]); + EXPECT_FLOAT_EQ(outvec[1][1],outstvec[1][1]); // stan::math::vector_d x3(3); // x3 << -1.0, 1.0, 10.0; diff --git a/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp b/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp index 8dbb32003c1..4ee758f5f0d 100644 --- a/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp +++ b/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp @@ -14,13 +14,32 @@ void test_log_sum_exp(const Eigen::Matrix& as) { } TEST(MathFunctions, log_sum_exp_mat) { + using stan::math::log_sum_exp; using Eigen::Dynamic; using Eigen::Matrix; - using stan::math::log_sum_exp; - Matrix m(3, 2); - m << 1, 2, 3, 4, 5, 6; - test_log_sum_exp(m); + Matrix m1(3, 2); + m1 << 1, 2, 3, 4, 5, 6; + std::vector st1{1, 2, 3, 4, 5, 6}; + Matrix m2(3, 2); + m2 << -1, -2, -3, -4, -5, -6; + std::vector st2{-1, -2, -3, -4, -5, -6}; + double m1_out = log_sum_exp(m1); + double m2_out = log_sum_exp(m2); + double st1_out = log_sum_exp(st1); + double st2_out = log_sum_exp(st2); + EXPECT_FLOAT_EQ(m1_out, st1_out); + EXPECT_FLOAT_EQ(m2_out, st2_out); + //test_log_sum_exp(m); + + std::vector st_m{m1,m2}; + std::vector> st_st{st1,st2}; + std::vector m_out = log_sum_exp(st_m); + std::vector st_out = log_sum_exp(st_st); + EXPECT_FLOAT_EQ(m1_out, m_out[0]); + EXPECT_FLOAT_EQ(m2_out, m_out[1]); + EXPECT_FLOAT_EQ(m1_out, st_out[0]); + EXPECT_FLOAT_EQ(m2_out, st_out[1]); Matrix v(3); v << 1, 2, 3; From 6b54c52956d7c49a18d154fcf2f819d2aa6456e0 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Mon, 16 Dec 2019 17:36:43 +0800 Subject: [PATCH 02/22] Add forwarding, rev & fwd versions --- stan/math/fwd/mat/fun/log_softmax.hpp | 60 +++++++++---------- stan/math/prim/mat/fun/head.hpp | 6 +- stan/math/prim/mat/fun/log_softmax.hpp | 8 +-- stan/math/prim/mat/fun/log_sum_exp.hpp | 6 +- .../prim/mat/vectorize/apply_vector_unary.hpp | 4 +- stan/math/rev/core/matrix_vari.hpp | 10 ++-- stan/math/rev/mat/fun/log_softmax.hpp | 58 +++++++++--------- stan/math/rev/mat/fun/log_sum_exp.hpp | 30 +++------- 8 files changed, 84 insertions(+), 98 deletions(-) diff --git a/stan/math/fwd/mat/fun/log_softmax.hpp b/stan/math/fwd/mat/fun/log_softmax.hpp index 711fa220b4e..e88f0ffaba0 100644 --- a/stan/math/fwd/mat/fun/log_softmax.hpp +++ b/stan/math/fwd/mat/fun/log_softmax.hpp @@ -6,44 +6,42 @@ #include #include #include +#include namespace stan { namespace math { -template -inline Eigen::Matrix, Eigen::Dynamic, 1> log_softmax( - const Eigen::Matrix, Eigen::Dynamic, 1>& alpha) { - using Eigen::Dynamic; - using Eigen::Matrix; - - Matrix alpha_t(alpha.size()); - for (int k = 0; k < alpha.size(); ++k) { - alpha_t(k) = alpha(k).val_; - } - - Matrix softmax_alpha_t = softmax(alpha_t); - Matrix log_softmax_alpha_t = log_softmax(alpha_t); - - Matrix, Dynamic, 1> log_softmax_alpha(alpha.size()); - for (int k = 0; k < alpha.size(); ++k) { - log_softmax_alpha(k).val_ = log_softmax_alpha_t(k); - log_softmax_alpha(k).d_ = 0; - } - - for (int m = 0; m < alpha.size(); ++m) { - T negative_alpha_m_d_times_softmax_alpha_t_m - = -alpha(m).d_ * softmax_alpha_t(m); - for (int k = 0; k < alpha.size(); ++k) { - if (m == k) { - log_softmax_alpha(k).d_ - += alpha(m).d_ + negative_alpha_m_d_times_softmax_alpha_t_m; - } else { - log_softmax_alpha(k).d_ += negative_alpha_m_d_times_softmax_alpha_t_m; +template >>...> +inline auto log_softmax(T&& x) { + return apply_vector_unary::apply(std::forward(x), [](auto& alpha){ + using Eigen::Dynamic; + using Eigen::Matrix; + + using T_scalar = value_type_t; + using T_fvar = typename T_scalar::Scalar; + + Matrix alpha_t = alpha.val(); + Matrix softmax_alpha_t = softmax(alpha_t); + + Matrix log_softmax_alpha(alpha.size()); + log_softmax_alpha.val() = log_softmax(alpha_t); + log_softmax_alpha.d().setZero(); + + for (int m = 0; m < alpha.size(); ++m) { + T_fvar negative_alpha_m_d_times_softmax_alpha_t_m + = -alpha(m).d_ * softmax_alpha_t(m); + for (int k = 0; k < alpha.size(); ++k) { + if (m == k) { + log_softmax_alpha(k).d_ + += alpha(m).d_ + negative_alpha_m_d_times_softmax_alpha_t_m; + } else { + log_softmax_alpha(k).d_ += negative_alpha_m_d_times_softmax_alpha_t_m; + } } } - } - return log_softmax_alpha; + return log_softmax_alpha; + }); } } // namespace math diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index 9d077bdd89e..456d11518fc 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -31,8 +31,8 @@ inline Eigen::Matrix head( }*/ template -inline auto head(const T& x, const T2& n) { - return apply_vector_unary::apply_scalar(x, n, [](auto& v, auto& m){ +inline auto head(T&& x, const T2& n) { + return apply_vector_unary::apply_scalar(std::forward(x), n, [](auto& v, auto& m){ if (m != 0){ if (v.rows() == 1){ check_column_index("head", "n", v, m); @@ -40,7 +40,7 @@ inline auto head(const T& x, const T2& n) { check_row_index("head", "n", v, m); } } - return v.head(m); + return v.head(m).eval(); }); } diff --git a/stan/math/prim/mat/fun/log_softmax.hpp b/stan/math/prim/mat/fun/log_softmax.hpp index 9dfc390ed71..5541f18c19a 100644 --- a/stan/math/prim/mat/fun/log_softmax.hpp +++ b/stan/math/prim/mat/fun/log_softmax.hpp @@ -38,11 +38,11 @@ namespace math { * @param[in] v Vector to transform. * @return Unit simplex result of the softmax transform of the vector. */ -template -inline auto log_softmax(const T& x) { - return apply_vector_unary::apply(x, [](auto& v){ +template >>...> +inline auto log_softmax(T&& x) { + return apply_vector_unary::apply(std::forward(x), [](auto& v){ check_nonzero_size("log_softmax", "v", v); - return v.array() - log_sum_exp(v); + return (v.array() - log_sum_exp(v)).matrix().eval(); }); } } // namespace math diff --git a/stan/math/prim/mat/fun/log_sum_exp.hpp b/stan/math/prim/mat/fun/log_sum_exp.hpp index 07ed5ad1910..56352699934 100644 --- a/stan/math/prim/mat/fun/log_sum_exp.hpp +++ b/stan/math/prim/mat/fun/log_sum_exp.hpp @@ -38,9 +38,9 @@ double log_sum_exp(const Eigen::Matrix& x) { return max + std::log((x.array() - max).exp().sum()); }*/ -template -inline auto log_sum_exp(const T& x) { - return apply_vector_unary::reduce(x, [](auto& v){ +template >>...> +inline auto log_sum_exp(T&& x) { + return apply_vector_unary::reduce(std::forward(x), [](auto& v){ if (v.size() == 0) { return -std::numeric_limits::infinity(); } diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp index df6c5e17c1b..d5a4b68f8de 100644 --- a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -30,7 +30,7 @@ struct apply_vector_unary> { }; template -struct apply_vector_unary> { +struct apply_vector_unary> { using scalar_t = scalar_type_t; using map_t = typename Eigen::Map>; @@ -59,7 +59,7 @@ struct apply_vector_unary> { template struct apply_vector_unary> { - using eigen_t = typename T::value_type; + using eigen_t = value_type_t; using scalar_t = typename eigen_t::Scalar; using return_t = std::vector #include +#include #include #include +#include namespace stan { namespace math { @@ -15,13 +17,11 @@ class op_matrix_vari : public vari { vari** vis_; public: - template - op_matrix_vari(double f, const Eigen::Matrix& vs) + template >>...> + op_matrix_vari(double f, T&& vs) : vari(f), size_(vs.size()) { vis_ = reinterpret_cast(operator new(sizeof(vari*) * vs.size())); - for (int i = 0; i < vs.size(); ++i) { - vis_[i] = vs(i).vi_; - } + Eigen::Map(vis_, vs.rows(), vs.cols()) = vs.vi(); } vari* operator[](size_t n) const { return vis_[n]; } size_t size() { return size_; } diff --git a/stan/math/rev/mat/fun/log_softmax.hpp b/stan/math/rev/mat/fun/log_softmax.hpp index ec990c1d673..f2bf6a2b6e8 100644 --- a/stan/math/rev/mat/fun/log_softmax.hpp +++ b/stan/math/rev/mat/fun/log_softmax.hpp @@ -54,41 +54,43 @@ class log_softmax_elt_vari : public vari { * @return Softmax of the input. * @throw std::domain_error If the input vector is size 0. */ -inline Eigen::Matrix log_softmax( - const Eigen::Matrix& alpha) { - const int a_size = alpha.size(); +template >>...> +inline auto log_softmax(T&& x) { + return apply_vector_unary::apply(std::forward(x), [](auto& alpha){ + const int a_size = alpha.size(); - check_nonzero_size("log_softmax", "alpha", alpha); + check_nonzero_size("log_softmax", "alpha", alpha); - // TODO(carpenter): replace with array alloc - vari** alpha_vi_array - = reinterpret_cast(vari::operator new(sizeof(vari*) * a_size)); - Eigen::Map(alpha_vi_array, a_size) = alpha.vi(); + // TODO(carpenter): replace with array alloc + vari** alpha_vi_array + = reinterpret_cast(vari::operator new(sizeof(vari*) * a_size)); + Eigen::Map(alpha_vi_array, a_size) = alpha.vi(); - vector_d alpha_d = alpha.val(); + vector_d alpha_d = alpha.val(); - // fold logic of math::softmax() and math::log_softmax() - // to save computations + // fold logic of math::softmax() and math::log_softmax() + // to save computations - vector_d diff = (alpha_d.array() - alpha_d.maxCoeff()); - vector_d softmax_alpha_d = diff.array().exp(); - double sum = softmax_alpha_d.sum(); - softmax_alpha_d.array() /= sum; - vector_d log_softmax_alpha_d = diff.array() - std::log(sum); + vector_d diff = (alpha_d.array() - alpha_d.maxCoeff()); + vector_d softmax_alpha_d = diff.array().exp(); + double sum = softmax_alpha_d.sum(); + softmax_alpha_d.array() /= sum; + vector_d log_softmax_alpha_d = diff.array() - std::log(sum); - // end fold - // TODO(carpenter): replace with array alloc - double* softmax_alpha_d_array - = reinterpret_cast(vari::operator new(sizeof(double) * a_size)); - Eigen::Map(softmax_alpha_d_array, a_size) = softmax_alpha_d; + // end fold + // TODO(carpenter): replace with array alloc + double* softmax_alpha_d_array + = reinterpret_cast(vari::operator new(sizeof(double) * a_size)); + Eigen::Map(softmax_alpha_d_array, a_size) = softmax_alpha_d; - vector_v log_softmax_alpha(a_size); - for (int k = 0; k < a_size; ++k) { - log_softmax_alpha(k) = var(new internal::log_softmax_elt_vari( - log_softmax_alpha_d[k], alpha_vi_array, softmax_alpha_d_array, a_size, - k)); - } - return log_softmax_alpha; + vector_v log_softmax_alpha(a_size); + for (int k = 0; k < a_size; ++k) { + log_softmax_alpha(k) = var(new internal::log_softmax_elt_vari( + log_softmax_alpha_d[k], alpha_vi_array, softmax_alpha_d_array, a_size, + k)); + } + return log_softmax_alpha; + }); } } // namespace math diff --git a/stan/math/rev/mat/fun/log_sum_exp.hpp b/stan/math/rev/mat/fun/log_sum_exp.hpp index 1e537f66c8c..2fcbe696db3 100644 --- a/stan/math/rev/mat/fun/log_sum_exp.hpp +++ b/stan/math/rev/mat/fun/log_sum_exp.hpp @@ -14,27 +14,11 @@ namespace math { namespace internal { -// these function and the following class just translate -// log_sum_exp for std::vector for Eigen::Matrix - -template -inline double log_sum_exp_as_double(const Eigen::Matrix& x) { - if (x.size() == 0) { - return -std::numeric_limits::infinity(); - } - - const double max = x.val().maxCoeff(); - if (!std::isfinite(max)) { - return max; - } - return max + std::log((x.val().array() - max).exp().sum()); -} - class log_sum_exp_matrix_vari : public op_matrix_vari { public: - template - explicit log_sum_exp_matrix_vari(const Eigen::Matrix& x) - : op_matrix_vari(log_sum_exp_as_double(x), x) {} + template + explicit log_sum_exp_matrix_vari(T&& x) + : op_matrix_vari(log_sum_exp(std::forward(x).val()), std::forward(x)) {} void chain() { Eigen::Map vis_map(vis_, size_); vis_map.adj().array() += adj_ * (vis_map.val().array() - val_).exp(); @@ -47,9 +31,11 @@ class log_sum_exp_matrix_vari : public op_matrix_vari { * * @param x matrix */ -template -inline var log_sum_exp(const Eigen::Matrix& x) { - return var(new internal::log_sum_exp_matrix_vari(x)); +template >>...> +inline auto log_sum_exp(T&& x) { + return apply_vector_unary::reduce(std::forward(x), [](auto& v){ + return var(new internal::log_sum_exp_matrix_vari(v)); + }); } } // namespace math From 2b085d151db69663bf3e197f59543fada6a7d6f5 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 17 Dec 2019 10:11:24 +0800 Subject: [PATCH 03/22] Add autodiff tests, remove arr versions --- stan/math/fwd/arr.hpp | 1 - stan/math/fwd/arr/fun/log_sum_exp.hpp | 31 ----------- stan/math/fwd/mat/fun/log_sum_exp.hpp | 17 ++++-- stan/math/prim/arr.hpp | 1 - stan/math/prim/arr/fun/log_sum_exp.hpp | 50 ----------------- stan/math/prim/mat/fun/head.hpp | 53 +----------------- stan/math/prim/mat/fun/log_sum_exp.hpp | 14 ----- .../prim/mat/prob/categorical_logit_lpmf.hpp | 2 +- stan/math/rev/arr.hpp | 1 - stan/math/rev/arr/fun/log_sum_exp.hpp | 55 ------------------- .../math/mix/mat/fun/log_softmax_test.cpp | 41 ++++++++++++++ .../math/mix/mat/fun/log_sum_exp_test.cpp | 6 +- 12 files changed, 57 insertions(+), 215 deletions(-) delete mode 100644 stan/math/fwd/arr/fun/log_sum_exp.hpp delete mode 100644 stan/math/prim/arr/fun/log_sum_exp.hpp delete mode 100644 stan/math/rev/arr/fun/log_sum_exp.hpp diff --git a/stan/math/fwd/arr.hpp b/stan/math/fwd/arr.hpp index 9cb494186b2..095f46a5ff5 100644 --- a/stan/math/fwd/arr.hpp +++ b/stan/math/fwd/arr.hpp @@ -7,7 +7,6 @@ #include #include -#include #include #include diff --git a/stan/math/fwd/arr/fun/log_sum_exp.hpp b/stan/math/fwd/arr/fun/log_sum_exp.hpp deleted file mode 100644 index 2e8fc3ddf83..00000000000 --- a/stan/math/fwd/arr/fun/log_sum_exp.hpp +++ /dev/null @@ -1,31 +0,0 @@ -#ifndef STAN_MATH_FWD_ARR_FUN_LOG_SUM_EXP_HPP -#define STAN_MATH_FWD_ARR_FUN_LOG_SUM_EXP_HPP - -#include -#include -#include -#include - -namespace stan { -namespace math { - -template -fvar log_sum_exp(const std::vector >& v) { - using std::exp; - std::vector vals(v.size()); - for (size_t i = 0; i < v.size(); ++i) { - vals[i] = v[i].val_; - } - T deriv(0.0); - T denominator(0.0); - for (size_t i = 0; i < v.size(); ++i) { - T exp_vi = exp(vals[i]); - denominator += exp_vi; - deriv += v[i].d_ * exp_vi; - } - return fvar(log_sum_exp(vals), deriv / denominator); -} - -} // namespace math -} // namespace stan -#endif diff --git a/stan/math/fwd/mat/fun/log_sum_exp.hpp b/stan/math/fwd/mat/fun/log_sum_exp.hpp index 8630737f4d9..4b995f1858c 100644 --- a/stan/math/fwd/mat/fun/log_sum_exp.hpp +++ b/stan/math/fwd/mat/fun/log_sum_exp.hpp @@ -8,13 +8,18 @@ namespace stan { namespace math { -template -fvar log_sum_exp(const Eigen::Matrix, R, C>& v) { - Eigen::Matrix vals = v.val(); - Eigen::Matrix exp_vals = vals.array().exp(); +template >>...> +inline auto log_sum_exp(T&& x) { + return apply_vector_unary::reduce(std::forward(x), [](auto& v){ + using T_scalar = value_type_t; + using T_fvar = typename T_scalar::Scalar; + using mat_type = Eigen::Matrix; + mat_type vals = v.val(); + mat_type exp_vals = vals.array().exp(); - return fvar(log_sum_exp(vals), - v.d().cwiseProduct(exp_vals).sum() / exp_vals.sum()); + return fvar(log_sum_exp(vals), + v.d().cwiseProduct(exp_vals).sum() / exp_vals.sum()); + }); } } // namespace math diff --git a/stan/math/prim/arr.hpp b/stan/math/prim/arr.hpp index 005c0d97b6f..6b2c0394d62 100644 --- a/stan/math/prim/arr.hpp +++ b/stan/math/prim/arr.hpp @@ -15,7 +15,6 @@ #include #include #include -#include #include #include #include diff --git a/stan/math/prim/arr/fun/log_sum_exp.hpp b/stan/math/prim/arr/fun/log_sum_exp.hpp deleted file mode 100644 index 1b7f7f2ddbe..00000000000 --- a/stan/math/prim/arr/fun/log_sum_exp.hpp +++ /dev/null @@ -1,50 +0,0 @@ -#ifndef STAN_MATH_PRIM_ARR_FUN_LOG_SUM_EXP_HPP -#define STAN_MATH_PRIM_ARR_FUN_LOG_SUM_EXP_HPP - -#include -#include -#include -#include -#include - -namespace stan { -namespace math { - -/** - * Return the log of the sum of the exponentiated values of the specified - * sequence of values. - * - * The function is defined as follows to prevent overflow in exponential - * calculations. - * - * \f$\log \sum_{n=1}^N \exp(x_n) = \max(x) + \log \sum_{n=1}^N \exp(x_n - - * \max(x))\f$. - * - * @param[in] x array of specified values - * @return The log of the sum of the exponentiated vector values. - */ -inline double log_sum_exp(const std::vector& x) { - using std::exp; - using std::log; - using std::numeric_limits; - double max = -numeric_limits::infinity(); - for (double xx : x) { - if (xx > max) { - max = xx; - } - } - - double sum = 0.0; - for (size_t ii = 0; ii < x.size(); ii++) { - if (x[ii] != -numeric_limits::infinity()) { - sum += exp(x[ii] - max); - } - } - - return max + log(sum); -} - -} // namespace math -} // namespace stan - -#endif diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index 456d11518fc..0ba3084839a 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -20,16 +20,7 @@ namespace math { * @param n Size of return. * @return The first n elements of v. * @throw std::out_of_range if n is out of range. - *//* -template -inline Eigen::Matrix head( - const Eigen::Matrix& v, size_t n) { - if (n != 0) { - check_row_index("head", "n", v, n); - } - return v.head(n); -}*/ - + */ template inline auto head(T&& x, const T2& n) { return apply_vector_unary::apply_scalar(std::forward(x), n, [](auto& v, auto& m){ @@ -44,48 +35,6 @@ inline auto head(T&& x, const T2& n) { }); } -/** - * Return the specified number of elements as a row vector - * from the front of the specified row vector. - * - * @tparam T Type of value in vector. - * @param rv Row vector. - * @param n Size of return row vector. - * @return The first n elements of rv. - * @throw std::out_of_range if n is out of range. - *//* -template -inline Eigen::Matrix head( - const Eigen::Matrix& rv, size_t n) { - if (n != 0) { - check_column_index("head", "n", rv, n); - } - return rv.head(n); -}*/ - -/** - * Return the specified number of elements as a standard vector - * from the front of the specified standard vector. - * - * @tparam T Type of value in vector. - * @param sv Standard vector. - * @param n Size of return. - * @return The first n elements of sv. - * @throw std::out_of_range if n is out of range. - *//* -template -std::vector head(const std::vector& sv, size_t n) { - if (n != 0) { - check_std_vector_index("head", "n", sv, n); - } - - std::vector s; - for (size_t i = 0; i < n; ++i) { - s.push_back(sv[i]); - } - return s; -}*/ - } // namespace math } // namespace stan #endif diff --git a/stan/math/prim/mat/fun/log_sum_exp.hpp b/stan/math/prim/mat/fun/log_sum_exp.hpp index 56352699934..a9f0578124b 100644 --- a/stan/math/prim/mat/fun/log_sum_exp.hpp +++ b/stan/math/prim/mat/fun/log_sum_exp.hpp @@ -24,20 +24,6 @@ namespace math { * @param[in] x Matrix of specified values * @return The log of the sum of the exponentiated vector values. */ - /* -template -double log_sum_exp(const Eigen::Matrix& x) { - if (x.size() == 0) { - return -std::numeric_limits::infinity(); - } - - const double max = x.maxCoeff(); - if (!std::isfinite(max)) { - return max; - } - return max + std::log((x.array() - max).exp().sum()); -}*/ - template >>...> inline auto log_sum_exp(T&& x) { return apply_vector_unary::reduce(std::forward(x), [](auto& v){ diff --git a/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp b/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp index 90b049a4f36..b934a1d7e56 100644 --- a/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp +++ b/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp @@ -4,7 +4,7 @@ #include #include #include -#include +#include #include #include #include diff --git a/stan/math/rev/arr.hpp b/stan/math/rev/arr.hpp index 74eeefd9822..2bf0107c633 100644 --- a/stan/math/rev/arr.hpp +++ b/stan/math/rev/arr.hpp @@ -7,7 +7,6 @@ #include #include -#include #include #include #include diff --git a/stan/math/rev/arr/fun/log_sum_exp.hpp b/stan/math/rev/arr/fun/log_sum_exp.hpp deleted file mode 100644 index 3bbaf1be1ad..00000000000 --- a/stan/math/rev/arr/fun/log_sum_exp.hpp +++ /dev/null @@ -1,55 +0,0 @@ -#ifndef STAN_MATH_REV_ARR_FUN_LOG_SUM_EXP_HPP -#define STAN_MATH_REV_ARR_FUN_LOG_SUM_EXP_HPP - -#include -#include -#include -#include -#include -#include - -namespace stan { -namespace math { - -namespace internal { -inline double log_sum_exp_as_double(const std::vector& x) { - using std::exp; - using std::log; - using std::numeric_limits; - double max = -numeric_limits::infinity(); - for (size_t i = 0; i < x.size(); ++i) { - if (x[i] > max) { - max = x[i].val(); - } - } - double sum = 0.0; - for (size_t i = 0; i < x.size(); ++i) { - if (x[i] != -numeric_limits::infinity()) { - sum += exp(x[i].val() - max); - } - } - return max + log(sum); -} - -class log_sum_exp_vector_vari : public op_vector_vari { - public: - explicit log_sum_exp_vector_vari(const std::vector& x) - : op_vector_vari(log_sum_exp_as_double(x), x) {} - void chain() { - for (size_t i = 0; i < size_; ++i) { - vis_[i]->adj_ += adj_ * calculate_chain(vis_[i]->val_, val_); - } - } -}; -} // namespace internal - -/** - * Returns the log sum of exponentials. - */ -inline var log_sum_exp(const std::vector& x) { - return var(new internal::log_sum_exp_vector_vari(x)); -} - -} // namespace math -} // namespace stan -#endif diff --git a/test/unit/math/mix/mat/fun/log_softmax_test.cpp b/test/unit/math/mix/mat/fun/log_softmax_test.cpp index bdd4cdfd085..c6c7d457f10 100644 --- a/test/unit/math/mix/mat/fun/log_softmax_test.cpp +++ b/test/unit/math/mix/mat/fun/log_softmax_test.cpp @@ -24,4 +24,45 @@ TEST(MathMixMatFun, logSoftmax) { Eigen::VectorXd x3c(3); x3c << 2, 1, 1; stan::test::expect_ad(f, x3c); + + Eigen::RowVectorXd rx0(0); // error case + stan::test::expect_ad(f, rx0); + + Eigen::RowVectorXd rx1(1); + rx1 << 0; + stan::test::expect_ad(f, rx1); + + Eigen::RowVectorXd rx2(2); + rx2 << -1, 1; + stan::test::expect_ad(f, rx2); + + Eigen::RowVectorXd rx3(3); + rx3 << -1, 1, 10; + stan::test::expect_ad(f, rx3); + + Eigen::RowVectorXd rx3b(3); + rx3b << 0, 1, 2; + stan::test::expect_ad(f, rx3b); + + Eigen::RowVectorXd rx3c(3); + rx3c << 2, 1, 1; + stan::test::expect_ad(f, rx3c); + + std::vector stx0(0); // error case + stan::test::expect_ad(f, stx0); + + std::vector stx1{0}; + stan::test::expect_ad(f, stx1); + + std::vector stx2{-1, 1}; + stan::test::expect_ad(f, stx2); + + std::vector stx3{-1, 1, 10}; + stan::test::expect_ad(f, stx3); + + std::vector stx3b{0, 1, 2}; + stan::test::expect_ad(f, stx3b); + + std::vector stx3c{2, 1, 1}; + stan::test::expect_ad(f, stx3c); } diff --git a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp index 2933019c1cb..ee695c76967 100644 --- a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp +++ b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp @@ -39,11 +39,11 @@ TEST(MathMixMatFun, logSumExp) { stan::test::expect_ad(f, x23); std::vector a1{0}; - stan::test::expect_ad(f, a1); + stan::test::expect_ad(tols, f, a1); std::vector a2{5, 2}; - stan::test::expect_ad(f, a2); + stan::test::expect_ad(tols, f, a2); std::vector a4{1, 2, 3, 4}; - stan::test::expect_ad(f, a4); + stan::test::expect_ad(tols, f, a4); } From ac3d5e36ec5ff05813351202b81301780cd0f027 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 24 Dec 2019 11:17:55 +0800 Subject: [PATCH 04/22] Nested testing --- stan/math/fwd/mat/fun/log_softmax.hpp | 14 ++--- stan/math/rev/mat/fun/log_softmax.hpp | 12 ++-- .../math/rev/mat/fun/log_softmax_test.cpp | 56 +++++++++++++++++++ 3 files changed, 68 insertions(+), 14 deletions(-) create mode 100644 test/unit/math/rev/mat/fun/log_softmax_test.cpp diff --git a/stan/math/fwd/mat/fun/log_softmax.hpp b/stan/math/fwd/mat/fun/log_softmax.hpp index e88f0ffaba0..3fbfb974895 100644 --- a/stan/math/fwd/mat/fun/log_softmax.hpp +++ b/stan/math/fwd/mat/fun/log_softmax.hpp @@ -13,22 +13,22 @@ namespace math { template >>...> inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [](auto& alpha){ + return apply_vector_unary::apply(std::forward(x), [&](auto& alpha){ using Eigen::Dynamic; using Eigen::Matrix; - using T_scalar = value_type_t; - using T_fvar = typename T_scalar::Scalar; + using T_value_type = value_type_t; + using T_scalar_type = typename T_value_type::Scalar; - Matrix alpha_t = alpha.val(); - Matrix softmax_alpha_t = softmax(alpha_t); + Matrix alpha_t = alpha.val(); + Matrix softmax_alpha_t = softmax(alpha_t); - Matrix log_softmax_alpha(alpha.size()); + Matrix log_softmax_alpha(alpha.size()); log_softmax_alpha.val() = log_softmax(alpha_t); log_softmax_alpha.d().setZero(); for (int m = 0; m < alpha.size(); ++m) { - T_fvar negative_alpha_m_d_times_softmax_alpha_t_m + T_scalar_type negative_alpha_m_d_times_softmax_alpha_t_m = -alpha(m).d_ * softmax_alpha_t(m); for (int k = 0; k < alpha.size(); ++k) { if (m == k) { diff --git a/stan/math/rev/mat/fun/log_softmax.hpp b/stan/math/rev/mat/fun/log_softmax.hpp index f2bf6a2b6e8..46c40348ba1 100644 --- a/stan/math/rev/mat/fun/log_softmax.hpp +++ b/stan/math/rev/mat/fun/log_softmax.hpp @@ -56,14 +56,13 @@ class log_softmax_elt_vari : public vari { */ template >>...> inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [](auto& alpha){ + return apply_vector_unary::apply(std::forward(x), [&](auto& alpha){ const int a_size = alpha.size(); check_nonzero_size("log_softmax", "alpha", alpha); - // TODO(carpenter): replace with array alloc vari** alpha_vi_array - = reinterpret_cast(vari::operator new(sizeof(vari*) * a_size)); + = ChainableStack::instance_->memalloc_.alloc_array(a_size); Eigen::Map(alpha_vi_array, a_size) = alpha.vi(); vector_d alpha_d = alpha.val(); @@ -74,14 +73,13 @@ inline auto log_softmax(T&& x) { vector_d diff = (alpha_d.array() - alpha_d.maxCoeff()); vector_d softmax_alpha_d = diff.array().exp(); double sum = softmax_alpha_d.sum(); - softmax_alpha_d.array() /= sum; vector_d log_softmax_alpha_d = diff.array() - std::log(sum); // end fold - // TODO(carpenter): replace with array alloc double* softmax_alpha_d_array - = reinterpret_cast(vari::operator new(sizeof(double) * a_size)); - Eigen::Map(softmax_alpha_d_array, a_size) = softmax_alpha_d; + = ChainableStack::instance_->memalloc_.alloc_array(a_size); + Eigen::Map(softmax_alpha_d_array, a_size) + = softmax_alpha_d.array() / sum; vector_v log_softmax_alpha(a_size); for (int k = 0; k < a_size; ++k) { diff --git a/test/unit/math/rev/mat/fun/log_softmax_test.cpp b/test/unit/math/rev/mat/fun/log_softmax_test.cpp new file mode 100644 index 00000000000..0161558cd24 --- /dev/null +++ b/test/unit/math/rev/mat/fun/log_softmax_test.cpp @@ -0,0 +1,56 @@ +#include +#include +#include +#include + +TEST(AgradRevMatrix, log_softmax_nested) { + using stan::math::log_softmax; + using stan::math::row_vector_v; + using stan::math::vector_v; + using stan::math::vector_d; + using stan::math::var; + + std::vector v(5); + std::vector rv(5); + std::vector> stv(5); + + std::vector v_result_iter(5); + std::vector rv_result_iter(5); + std::vector> stv_result_iter(5); + + for(int i = 0; i < 5; ++i){ + v[i] = vector_d::Random(10); + v[i].adj() = vector_d::Random(10); + rv[i] = v[i].transpose(); + stv[i] = std::vector(v[i].data(), v[i].data() + v[i].size()); + + v_result_iter[i] = log_softmax(v[i]); + rv_result_iter[i] = log_softmax(rv[i]); + stv_result_iter[i] = log_softmax(stv[i]); + } + + std::vector v_result = log_softmax(v); + std::vector rv_result = log_softmax(rv); + std::vector> stv_result = log_softmax(stv); + + for(int i = 0; i < 5; ++i){ + Eigen::Map stv_result_map(stv_result[i].data(), + stv_result[i].size()); + Eigen::Map stv_result_iter_map(stv_result_iter[i].data(), + stv_result_iter[i].size()); + expect_matrix_eq(v_result[i].val().array() * 3, + v_result_iter[i].val() + + rv_result_iter[i].val().transpose() + + stv_result_iter_map.val()); + expect_matrix_eq(v_result[i].adj().array() * 3, + v_result_iter[i].adj() + + rv_result_iter[i].adj().transpose() + + stv_result_iter_map.adj()); + expect_matrix_eq(v_result[i].val().array() * 2, + rv_result[i].val().transpose() + + stv_result_map.val()); + expect_matrix_eq(v_result[i].adj().array() * 2, + rv_result[i].adj().transpose() + + stv_result_map.adj()); + } +} From 77f3a59dcd66808fa8d2cdb6c7ba6f972f2a0a00 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Thu, 26 Dec 2019 21:08:41 +0800 Subject: [PATCH 05/22] Fix tests, update doc --- stan/math/fwd/mat/fun/log_softmax.hpp | 25 ++++++--- stan/math/fwd/mat/fun/log_sum_exp.hpp | 25 +++++++-- stan/math/prim/mat/fun/log_softmax.hpp | 6 +- stan/math/prim/mat/fun/log_sum_exp.hpp | 7 ++- stan/math/rev/mat/fun/log_softmax.hpp | 4 +- stan/math/rev/mat/fun/log_sum_exp.hpp | 7 ++- .../math/mix/mat/fun/log_softmax_test.cpp | 38 ++++++++++--- .../math/mix/mat/fun/log_sum_exp_test.cpp | 21 ++++--- .../math/rev/mat/fun/log_softmax_test.cpp | 56 ------------------- 9 files changed, 95 insertions(+), 94 deletions(-) delete mode 100644 test/unit/math/rev/mat/fun/log_softmax_test.cpp diff --git a/stan/math/fwd/mat/fun/log_softmax.hpp b/stan/math/fwd/mat/fun/log_softmax.hpp index 3fbfb974895..ad6b1cca21b 100644 --- a/stan/math/fwd/mat/fun/log_softmax.hpp +++ b/stan/math/fwd/mat/fun/log_softmax.hpp @@ -11,24 +11,31 @@ namespace stan { namespace math { +/** + * Return the softmax of the specified vector or container of vectors. + * Softmax is guaranteed to return a simplex. + * + * @tparam T Type of input vector or matrix. + * @param[in] x Unconstrained input vector. + * @return Softmax of the input. + * @throw std::domain_error If the input vector is size 0. + */ template >>...> inline auto log_softmax(T&& x) { return apply_vector_unary::apply(std::forward(x), [&](auto& alpha){ - using Eigen::Dynamic; - using Eigen::Matrix; - using T_value_type = value_type_t; - using T_scalar_type = typename T_value_type::Scalar; + using T_fvar = value_type_t; + using T_fvar_inner = typename T_fvar::Scalar; + + Eigen::Matrix alpha_t = alpha.val(); + Eigen::Matrix softmax_alpha_t = softmax(alpha_t); - Matrix alpha_t = alpha.val(); - Matrix softmax_alpha_t = softmax(alpha_t); - - Matrix log_softmax_alpha(alpha.size()); + Eigen::Matrix log_softmax_alpha(alpha.size()); log_softmax_alpha.val() = log_softmax(alpha_t); log_softmax_alpha.d().setZero(); for (int m = 0; m < alpha.size(); ++m) { - T_scalar_type negative_alpha_m_d_times_softmax_alpha_t_m + T_fvar_inner negative_alpha_m_d_times_softmax_alpha_t_m = -alpha(m).d_ * softmax_alpha_t(m); for (int k = 0; k < alpha.size(); ++k) { if (m == k) { diff --git a/stan/math/fwd/mat/fun/log_sum_exp.hpp b/stan/math/fwd/mat/fun/log_sum_exp.hpp index 4b995f1858c..6eb936a9804 100644 --- a/stan/math/fwd/mat/fun/log_sum_exp.hpp +++ b/stan/math/fwd/mat/fun/log_sum_exp.hpp @@ -2,22 +2,37 @@ #define STAN_MATH_FWD_MAT_FUN_LOG_SUM_EXP_HPP #include +#include #include #include namespace stan { namespace math { +/** + * Return the log of the sum of the exponentiated values of the specified + * matrix of values. The matrix may be a full matrix, a vector, + * a row vector, or a container of these. + * + * The function is defined as follows to prevent overflow in exponential + * calculations. + * + * \f$\log \sum_{n=1}^N \exp(x_n) = \max(x) + \log \sum_{n=1}^N \exp(x_n - + * \max(x))\f$. + * + * @tparam T Type of input vector or matrix. + * @param[in] x Matrix of specified values. + * @return The log of the sum of the exponentiated vector values. + */ template >>...> inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [](auto& v){ - using T_scalar = value_type_t; - using T_fvar = typename T_scalar::Scalar; - using mat_type = Eigen::Matrix; + return apply_vector_unary::reduce(std::forward(x), [&](auto& v){ + using T_fvar_inner = typename value_type_t::Scalar; + using mat_type = Eigen::Matrix; mat_type vals = v.val(); mat_type exp_vals = vals.array().exp(); - return fvar(log_sum_exp(vals), + return fvar(log_sum_exp(vals), v.d().cwiseProduct(exp_vals).sum() / exp_vals.sum()); }); } diff --git a/stan/math/prim/mat/fun/log_softmax.hpp b/stan/math/prim/mat/fun/log_softmax.hpp index 5541f18c19a..91195d3ddd5 100644 --- a/stan/math/prim/mat/fun/log_softmax.hpp +++ b/stan/math/prim/mat/fun/log_softmax.hpp @@ -34,15 +34,15 @@ namespace math { * \right. * \f$ * - * @tparam T Scalar type of values in vector. + * @tparam T Type of input vector to transform. * @param[in] v Vector to transform. * @return Unit simplex result of the softmax transform of the vector. */ template >>...> inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [](auto& v){ + return apply_vector_unary::apply(std::forward(x), [&](auto& v){ check_nonzero_size("log_softmax", "v", v); - return (v.array() - log_sum_exp(v)).matrix().eval(); + return (v.array() - log_sum_exp(v)).matrix(); }); } } // namespace math diff --git a/stan/math/prim/mat/fun/log_sum_exp.hpp b/stan/math/prim/mat/fun/log_sum_exp.hpp index a9f0578124b..fe6d25ce013 100644 --- a/stan/math/prim/mat/fun/log_sum_exp.hpp +++ b/stan/math/prim/mat/fun/log_sum_exp.hpp @@ -13,7 +13,7 @@ namespace math { /** * Return the log of the sum of the exponentiated values of the specified * matrix of values. The matrix may be a full matrix, a vector, - * or a row vector. + * a row vector, or a container of these. * * The function is defined as follows to prevent overflow in exponential * calculations. @@ -21,12 +21,13 @@ namespace math { * \f$\log \sum_{n=1}^N \exp(x_n) = \max(x) + \log \sum_{n=1}^N \exp(x_n - * \max(x))\f$. * - * @param[in] x Matrix of specified values + * @tparam T Type of input vector or matrix. + * @param[in] x Matrix of specified values. * @return The log of the sum of the exponentiated vector values. */ template >>...> inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [](auto& v){ + return apply_vector_unary::reduce(std::forward(x), [&](auto& v){ if (v.size() == 0) { return -std::numeric_limits::infinity(); } diff --git a/stan/math/rev/mat/fun/log_softmax.hpp b/stan/math/rev/mat/fun/log_softmax.hpp index 46c40348ba1..64834de07ff 100644 --- a/stan/math/rev/mat/fun/log_softmax.hpp +++ b/stan/math/rev/mat/fun/log_softmax.hpp @@ -2,6 +2,7 @@ #define STAN_MATH_REV_MAT_FUN_LOG_SOFTMAX_HPP #include +#include #include #include #include @@ -50,7 +51,8 @@ class log_softmax_elt_vari : public vari { * * The gradient calculations are unfolded. * - * @param alpha Unconstrained input vector. + * @tparam T Type of input vector or matrix. + * @param[in] x Unconstrained input vector. * @return Softmax of the input. * @throw std::domain_error If the input vector is size 0. */ diff --git a/stan/math/rev/mat/fun/log_sum_exp.hpp b/stan/math/rev/mat/fun/log_sum_exp.hpp index 2fcbe696db3..fcaf8da39af 100644 --- a/stan/math/rev/mat/fun/log_sum_exp.hpp +++ b/stan/math/rev/mat/fun/log_sum_exp.hpp @@ -5,6 +5,7 @@ #include #include #include +#include #include #include #include @@ -18,7 +19,8 @@ class log_sum_exp_matrix_vari : public op_matrix_vari { public: template explicit log_sum_exp_matrix_vari(T&& x) - : op_matrix_vari(log_sum_exp(std::forward(x).val()), std::forward(x)) {} + : op_matrix_vari(log_sum_exp(std::forward(x).val()), + std::forward(x)) {} void chain() { Eigen::Map vis_map(vis_, size_); vis_map.adj().array() += adj_ * (vis_map.val().array() - val_).exp(); @@ -29,11 +31,12 @@ class log_sum_exp_matrix_vari : public op_matrix_vari { /** * Returns the log sum of exponentials. * + * @tparam T Type of input vector or matrix. * @param x matrix */ template >>...> inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [](auto& v){ + return apply_vector_unary::reduce(std::forward(x), [&](auto& v){ return var(new internal::log_sum_exp_matrix_vari(v)); }); } diff --git a/test/unit/math/mix/mat/fun/log_softmax_test.cpp b/test/unit/math/mix/mat/fun/log_softmax_test.cpp index c6c7d457f10..db6a51ac74a 100644 --- a/test/unit/math/mix/mat/fun/log_softmax_test.cpp +++ b/test/unit/math/mix/mat/fun/log_softmax_test.cpp @@ -1,7 +1,8 @@ #include TEST(MathMixMatFun, logSoftmax) { - auto f = [](const auto& x) { return stan::math::log_softmax(x); }; + auto f = [](const auto& x) { return stan::math::log_softmax(x).eval(); }; + //Column Vectors Eigen::VectorXd x0(0); // error case stan::test::expect_ad(f, x0); @@ -25,6 +26,7 @@ TEST(MathMixMatFun, logSoftmax) { x3c << 2, 1, 1; stan::test::expect_ad(f, x3c); + //Row Vectors Eigen::RowVectorXd rx0(0); // error case stan::test::expect_ad(f, rx0); @@ -48,21 +50,43 @@ TEST(MathMixMatFun, logSoftmax) { rx3c << 2, 1, 1; stan::test::expect_ad(f, rx3c); + auto g = [](const auto& x) { return stan::math::log_softmax(x); }; + + //std vectors std::vector stx0(0); // error case - stan::test::expect_ad(f, stx0); + stan::test::expect_ad(g, stx0); std::vector stx1{0}; - stan::test::expect_ad(f, stx1); + stan::test::expect_ad(g, stx1); std::vector stx2{-1, 1}; - stan::test::expect_ad(f, stx2); + stan::test::expect_ad(g, stx2); std::vector stx3{-1, 1, 10}; - stan::test::expect_ad(f, stx3); + stan::test::expect_ad(g, stx3); std::vector stx3b{0, 1, 2}; - stan::test::expect_ad(f, stx3b); + stan::test::expect_ad(g, stx3b); std::vector stx3c{2, 1, 1}; - stan::test::expect_ad(f, stx3c); + stan::test::expect_ad(g, stx3c); + + //Nested containers + std::vector stvx0{x0, x0}; // error case + stan::test::expect_ad(g, stvx0); + + std::vector stvx1{x1, x1}; + stan::test::expect_ad(g, stvx1); + + std::vector strx0{rx0, rx0}; // error case + stan::test::expect_ad(g, strx0); + + std::vector strx1{rx1, rx1}; + stan::test::expect_ad(g, strx1); + + std::vector> ststx0{stx0, stx0}; // error case + stan::test::expect_ad(g, ststx0); + + std::vector> ststx1{stx1, stx1}; + stan::test::expect_ad(g, ststx1); } diff --git a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp index ee695c76967..dbfb9b1204d 100644 --- a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp +++ b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp @@ -32,18 +32,23 @@ TEST(MathMixMatFun, logSumExp) { stan::test::expect_ad(tols, f, x); Eigen::RowVectorXd rx = x; stan::test::expect_ad(tols, f, rx); + std::vector stx = std::vector(x.data(), x.data() + x.size()); + stan::test::expect_ad(tols, f, stx); } + Eigen::MatrixXd x23(2, 2); x23 << 1, 2, 3, 4; stan::test::expect_ad(f, x23); - std::vector a1{0}; - stan::test::expect_ad(tols, f, a1); - - std::vector a2{5, 2}; - stan::test::expect_ad(tols, f, a2); - - std::vector a4{1, 2, 3, 4}; - stan::test::expect_ad(tols, f, a4); + std::vector stvx{x2, x2b, x2c}; + stan::test::expect_ad(tols, f, stvx); + std::vector strx{x2.transpose(), + x2b.transpose(), + x2c.transpose()}; + stan::test::expect_ad(tols, f, strx); + std::vector> ststx{std::vector(x2.data(), x2.data() + x2.size()), + std::vector(x2b.data(), x2b.data() + x2b.size()), + std::vector(x2c.data(), x2c.data() + x2c.size())}; + stan::test::expect_ad(tols, f, ststx); } diff --git a/test/unit/math/rev/mat/fun/log_softmax_test.cpp b/test/unit/math/rev/mat/fun/log_softmax_test.cpp deleted file mode 100644 index 0161558cd24..00000000000 --- a/test/unit/math/rev/mat/fun/log_softmax_test.cpp +++ /dev/null @@ -1,56 +0,0 @@ -#include -#include -#include -#include - -TEST(AgradRevMatrix, log_softmax_nested) { - using stan::math::log_softmax; - using stan::math::row_vector_v; - using stan::math::vector_v; - using stan::math::vector_d; - using stan::math::var; - - std::vector v(5); - std::vector rv(5); - std::vector> stv(5); - - std::vector v_result_iter(5); - std::vector rv_result_iter(5); - std::vector> stv_result_iter(5); - - for(int i = 0; i < 5; ++i){ - v[i] = vector_d::Random(10); - v[i].adj() = vector_d::Random(10); - rv[i] = v[i].transpose(); - stv[i] = std::vector(v[i].data(), v[i].data() + v[i].size()); - - v_result_iter[i] = log_softmax(v[i]); - rv_result_iter[i] = log_softmax(rv[i]); - stv_result_iter[i] = log_softmax(stv[i]); - } - - std::vector v_result = log_softmax(v); - std::vector rv_result = log_softmax(rv); - std::vector> stv_result = log_softmax(stv); - - for(int i = 0; i < 5; ++i){ - Eigen::Map stv_result_map(stv_result[i].data(), - stv_result[i].size()); - Eigen::Map stv_result_iter_map(stv_result_iter[i].data(), - stv_result_iter[i].size()); - expect_matrix_eq(v_result[i].val().array() * 3, - v_result_iter[i].val() + - rv_result_iter[i].val().transpose() + - stv_result_iter_map.val()); - expect_matrix_eq(v_result[i].adj().array() * 3, - v_result_iter[i].adj() + - rv_result_iter[i].adj().transpose() + - stv_result_iter_map.adj()); - expect_matrix_eq(v_result[i].val().array() * 2, - rv_result[i].val().transpose() + - stv_result_map.val()); - expect_matrix_eq(v_result[i].adj().array() * 2, - rv_result[i].adj().transpose() + - stv_result_map.adj()); - } -} From b3a11326fd93ad8266af8dbd35af3cc84c714173 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Fri, 27 Dec 2019 13:41:43 +0800 Subject: [PATCH 06/22] Tidy doc --- stan/math/prim/mat/fun/head.hpp | 10 +- .../prim/mat/vectorize/apply_vector_unary.hpp | 125 +++++++++++++++--- stan/math/rev/core/matrix_vari.hpp | 2 +- test/unit/math/mix/mat/fun/head_test.cpp | 12 +- 4 files changed, 126 insertions(+), 23 deletions(-) diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index 0ba3084839a..f969518b54e 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -15,15 +15,17 @@ namespace math { * Return the specified number of elements as a vector * from the front of the specified vector. * - * @tparam T Type of value in vector. - * @param v Vector input. + * @tparam T Type of input vector. + * @tparam T2 Type of size variable. + * @param x Vector input. * @param n Size of return. * @return The first n elements of v. * @throw std::out_of_range if n is out of range. */ template inline auto head(T&& x, const T2& n) { - return apply_vector_unary::apply_scalar(std::forward(x), n, [](auto& v, auto& m){ + return apply_vector_unary::apply_scalar(std::forward(x), n, + [&](auto& v, auto& m){ if (m != 0){ if (v.rows() == 1){ check_column_index("head", "n", v, m); @@ -31,7 +33,7 @@ inline auto head(T&& x, const T2& n) { check_row_index("head", "n", v, m); } } - return v.head(m).eval(); + return v.head(m); }); } diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp index d5a4b68f8de..10e35fdef17 100644 --- a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -8,49 +8,137 @@ namespace stan { namespace math { +//Forward declaration to allow specialisations template struct apply_vector_unary {}; +/** + * Base template class for vectorization of unary vector functions + * defined by applying a functor to a standard library vector, Eigen dense + * matrix expression template, or container of these. + * + * Three taxonomies of unary vector functions are implemented: + * - f(vector) -> vector + * - f(vector) -> scalar + * - f(vector, scalar) -> vector + * + * This base template class takes (and returns) Eigen expression templates. + */ template struct apply_vector_unary> { + /** + * Member function for applying a functor to a vector and subsequently + * returning a vector. The 'auto' return type is used here so that an + * expression template is returned. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x Eigen input to which operation is applied. + * @param f functor to apply to Eigen input. + * @return Eigen expression template with result of applying functor + * to input + */ template static inline auto apply(const T& x, const F& f) { return f(x); } + /** + * Member function for applying a functor to a vector and a scalar + * and subsequently returning a vector. The 'auto' return type is + * used here so that an expression template is returned. + * + * @tparam T Type of argument to which functor is applied. + * @tparam T2 Type of scalar to pass to functor. + * @tparam F Type of functor to apply. + * @param x Eigen input to which operation is applied. + * @param y scalar passed to functor. + * @param f functor to apply to Eigen input. + * @return Eigen expression template with result of applying functor + * and scalar to input + */ template static inline auto apply_scalar(const T& x, const T2& y, const F& f) { return f(x, y); } + /** + * Member function for applying a functor to a vector and subsequently + * returning a scalar. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x Eigen input to which operation is applied. + * @param f functor to apply to Eigen input. + * @return scalar result of applying functor to input. + */ template static inline auto reduce(const T& x, const F& f) { return f(x); } }; +/** + * Specialisation for use with (non-nested) std::vectors. Inputs are mapped + * to Eigen column vectors and then passed to the base (Eigen) template. + */ template struct apply_vector_unary> { using scalar_t = scalar_type_t; - using map_t = - typename Eigen::Map>; - + using map_t = typename Eigen::Map>; + + /** + * Member function for applying a functor to a vector and subsequently + * returning a vector. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x std::vector input to which operation is applied. + * @param f functor to apply to vector input. + * @return std::vector with result of applying functor to input. + */ template static inline std::vector apply(const T& x, const F& f) { - Eigen::Matrix result - = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); + Eigen::Matrix result + = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); return std::vector(result.data(), - result.data() + result.size()); + result.data() + result.size()); } + /** + * Member function for applying a functor to a vector and a scalar + * and subsequently returning a vector. The 'auto' return type is + * used here so that an expression template is returned. + * + * @tparam T Type of argument to which functor is applied. + * @tparam T2 Type of scalar to pass to functor. + * @tparam F Type of functor to apply. + * @param x std::vector input to which operation is applied. + * @param y scalar passed to functor. + * @param f functor to apply to vector input. + * @return std::vector with result of applying functor and scalar to input. + */ template - static inline std::vector apply_scalar(const T& x, const T2& y, const F& f) { - Eigen::Matrix result - = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), y, f); + static inline std::vector apply_scalar(const T& x, + const T2& y, + const F& f) { + Eigen::Matrix result + = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), + y, f); return std::vector(result.data(), result.data() + result.size()); } + /** + * Member function for applying a functor to a vector and subsequently + * returning a scalar. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x Eigen input to which operation is applied. + * @param f functor to apply to std::vector input. + * @return scalar result of applying functor to input vector. + */ template static inline scalar_t reduce(const T& x, const F& f) { return apply_vector_unary::reduce(as_column_vector_or_scalar(x), f); @@ -99,15 +187,16 @@ struct apply_vector_unary> { using scalar_t = scalar_type_t; using return_t = typename std::vector>; using map_t = - typename Eigen::Map>; + typename Eigen::Map>; template static inline return_t apply(const T& x, const F& f) { size_t x_size = x.size(); return_t result(x_size); - Eigen::Matrix inter; + Eigen::Matrix inter; for(size_t i = 0; i < x_size; ++i){ - inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), f); + inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), + f); result[i] = std::vector(inter.data(), inter.data() + inter.size()); } @@ -119,9 +208,10 @@ struct apply_vector_unary> { scalar_seq_view y_vec(y); size_t x_size = x.size(); return_t result(x_size); - Eigen::Matrix inter; + Eigen::Matrix inter; for(size_t i = 0; i < x_size; ++i){ - inter = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x[i]), y_vec[i], f); + inter = apply_vector_unary::apply_scalar( + as_column_vector_or_scalar(x[i]), y_vec[i], f); result[i] = std::vector(inter.data(), inter.data() + inter.size()); } @@ -132,9 +222,10 @@ struct apply_vector_unary> { static inline std::vector reduce(const T& x, const F& f) { size_t x_size = x.size(); std::vector result(x_size); - for(size_t i = 0; i < x_size; ++i) - result[i] = apply_vector_unary::reduce(as_column_vector_or_scalar(x[i]), f); - + for(size_t i = 0; i < x_size; ++i){ + result[i] = apply_vector_unary::reduce( + as_column_vector_or_scalar(x[i]), f); + } return result; } }; diff --git a/stan/math/rev/core/matrix_vari.hpp b/stan/math/rev/core/matrix_vari.hpp index 341c9f8de96..95f41e99c7d 100644 --- a/stan/math/rev/core/matrix_vari.hpp +++ b/stan/math/rev/core/matrix_vari.hpp @@ -20,7 +20,7 @@ class op_matrix_vari : public vari { template >>...> op_matrix_vari(double f, T&& vs) : vari(f), size_(vs.size()) { - vis_ = reinterpret_cast(operator new(sizeof(vari*) * vs.size())); + vis_ = ChainableStack::instance_->memalloc_.alloc_array(size_); Eigen::Map(vis_, vs.rows(), vs.cols()) = vs.vi(); } vari* operator[](size_t n) const { return vis_[n]; } diff --git a/test/unit/math/mix/mat/fun/head_test.cpp b/test/unit/math/mix/mat/fun/head_test.cpp index a5a4305d105..893a880ad87 100644 --- a/test/unit/math/mix/mat/fun/head_test.cpp +++ b/test/unit/math/mix/mat/fun/head_test.cpp @@ -2,16 +2,26 @@ #include auto f(int i) { + return [=](const auto& y) { return stan::math::head(y, i).eval(); }; +} + +auto g(int i) { return [=](const auto& y) { return stan::math::head(y, i); }; } template void expect_head(const T& x, int n) { + std::vector> stx{x, x, x}; Eigen::VectorXd v = stan::test::to_vector(x); + std::vector stv{v, v, v}; Eigen::RowVectorXd rv = stan::test::to_row_vector(x); - stan::test::expect_ad(f(n), x); + std::vector strv{rv, rv, rv}; + stan::test::expect_ad(g(n), x); + stan::test::expect_ad(g(n), stx); stan::test::expect_ad(f(n), v); + stan::test::expect_ad(g(n), stv); stan::test::expect_ad(f(n), rv); + stan::test::expect_ad(g(n), strv); } TEST(MathMixMatFun, head) { From a4b83a7dbda841f48fd4f19359964384285e4f6f Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Fri, 27 Dec 2019 23:10:33 +0800 Subject: [PATCH 07/22] Cpplint --- stan/math/fwd/mat/fun/log_softmax.hpp | 3 +-- stan/math/prim/mat/fun/head.hpp | 4 +-- .../prim/mat/prob/categorical_logit_lpmf.hpp | 1 - .../prim/mat/vectorize/apply_vector_unary.hpp | 27 ++++++++++--------- .../math/mix/mat/fun/log_softmax_test.cpp | 9 ++++--- .../math/mix/mat/fun/log_sum_exp_test.cpp | 11 +++++--- test/unit/math/prim/mat/fun/head_test.cpp | 2 +- .../math/prim/mat/fun/log_softmax_test.cpp | 21 ++++++++------- .../math/prim/mat/fun/log_sum_exp_test.cpp | 6 ++--- 9 files changed, 45 insertions(+), 39 deletions(-) diff --git a/stan/math/fwd/mat/fun/log_softmax.hpp b/stan/math/fwd/mat/fun/log_softmax.hpp index ad6b1cca21b..b40887776bf 100644 --- a/stan/math/fwd/mat/fun/log_softmax.hpp +++ b/stan/math/fwd/mat/fun/log_softmax.hpp @@ -23,10 +23,9 @@ namespace math { template >>...> inline auto log_softmax(T&& x) { return apply_vector_unary::apply(std::forward(x), [&](auto& alpha){ - using T_fvar = value_type_t; using T_fvar_inner = typename T_fvar::Scalar; - + Eigen::Matrix alpha_t = alpha.val(); Eigen::Matrix softmax_alpha_t = softmax(alpha_t); diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index f969518b54e..d32f429a5ec 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -26,8 +26,8 @@ template inline auto head(T&& x, const T2& n) { return apply_vector_unary::apply_scalar(std::forward(x), n, [&](auto& v, auto& m){ - if (m != 0){ - if (v.rows() == 1){ + if (m != 0) { + if (v.rows() == 1) { check_column_index("head", "n", v, m); } else { check_row_index("head", "n", v, m); diff --git a/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp b/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp index b934a1d7e56..409df3f7c77 100644 --- a/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp +++ b/stan/math/prim/mat/prob/categorical_logit_lpmf.hpp @@ -4,7 +4,6 @@ #include #include #include -#include #include #include #include diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp index 10e35fdef17..ee036196df0 100644 --- a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -8,7 +8,7 @@ namespace stan { namespace math { -//Forward declaration to allow specialisations +// Forward declaration to allow specialisations template struct apply_vector_unary {}; @@ -101,7 +101,7 @@ struct apply_vector_unary> { static inline std::vector apply(const T& x, const F& f) { Eigen::Matrix result = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); - return std::vector(result.data(), + return std::vector(result.data(), result.data() + result.size()); } @@ -125,8 +125,8 @@ struct apply_vector_unary> { Eigen::Matrix result = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), y, f); - return std::vector(result.data(), - result.data() + result.size()); + return std::vector(result.data(), + result.data() + result.size()); } /** @@ -145,6 +145,9 @@ struct apply_vector_unary> { } }; +/** + * Specialisation for use with std::vectors of Eigen types + */ template struct apply_vector_unary> { using eigen_t = value_type_t; @@ -157,7 +160,7 @@ struct apply_vector_unary> { static inline return_t apply(const T& x, const F& f) { size_t x_size = x.size(); return_t result(x_size); - for(size_t i = 0; i < x_size; ++i) + for (size_t i = 0; i < x_size; ++i) result[i] = apply_vector_unary::apply(x[i], f); return result; } @@ -167,7 +170,7 @@ struct apply_vector_unary> { scalar_seq_view y_vec(y); size_t x_size = x.size(); return_t result(x_size); - for(size_t i = 0; i < x_size; ++i) + for (size_t i = 0; i < x_size; ++i) result[i] = apply_vector_unary::apply_scalar(x[i], y_vec[i], f); return result; } @@ -176,7 +179,7 @@ struct apply_vector_unary> { static inline std::vector reduce(const T& x, const F& f) { size_t x_size = x.size(); std::vector result(x_size); - for(size_t i = 0; i < x_size; ++i) + for (size_t i = 0; i < x_size; ++i) result[i] = apply_vector_unary::reduce(x[i], f); return result; } @@ -194,10 +197,10 @@ struct apply_vector_unary> { size_t x_size = x.size(); return_t result(x_size); Eigen::Matrix inter; - for(size_t i = 0; i < x_size; ++i){ + for (size_t i = 0; i < x_size; ++i) { inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), f); - result[i] = std::vector(inter.data(), + result[i] = std::vector(inter.data(), inter.data() + inter.size()); } return result; @@ -209,10 +212,10 @@ struct apply_vector_unary> { size_t x_size = x.size(); return_t result(x_size); Eigen::Matrix inter; - for(size_t i = 0; i < x_size; ++i){ + for (size_t i = 0; i < x_size; ++i) { inter = apply_vector_unary::apply_scalar( as_column_vector_or_scalar(x[i]), y_vec[i], f); - result[i] = std::vector(inter.data(), + result[i] = std::vector(inter.data(), inter.data() + inter.size()); } return result; @@ -222,7 +225,7 @@ struct apply_vector_unary> { static inline std::vector reduce(const T& x, const F& f) { size_t x_size = x.size(); std::vector result(x_size); - for(size_t i = 0; i < x_size; ++i){ + for (size_t i = 0; i < x_size; ++i) { result[i] = apply_vector_unary::reduce( as_column_vector_or_scalar(x[i]), f); } diff --git a/test/unit/math/mix/mat/fun/log_softmax_test.cpp b/test/unit/math/mix/mat/fun/log_softmax_test.cpp index db6a51ac74a..bee6b4ed61e 100644 --- a/test/unit/math/mix/mat/fun/log_softmax_test.cpp +++ b/test/unit/math/mix/mat/fun/log_softmax_test.cpp @@ -1,8 +1,9 @@ #include +#include TEST(MathMixMatFun, logSoftmax) { auto f = [](const auto& x) { return stan::math::log_softmax(x).eval(); }; - //Column Vectors + // Column Vectors Eigen::VectorXd x0(0); // error case stan::test::expect_ad(f, x0); @@ -26,7 +27,7 @@ TEST(MathMixMatFun, logSoftmax) { x3c << 2, 1, 1; stan::test::expect_ad(f, x3c); - //Row Vectors + // Row Vectors Eigen::RowVectorXd rx0(0); // error case stan::test::expect_ad(f, rx0); @@ -52,7 +53,7 @@ TEST(MathMixMatFun, logSoftmax) { auto g = [](const auto& x) { return stan::math::log_softmax(x); }; - //std vectors + // std vectors std::vector stx0(0); // error case stan::test::expect_ad(g, stx0); @@ -71,7 +72,7 @@ TEST(MathMixMatFun, logSoftmax) { std::vector stx3c{2, 1, 1}; stan::test::expect_ad(g, stx3c); - //Nested containers + // Nested containers std::vector stvx0{x0, x0}; // error case stan::test::expect_ad(g, stvx0); diff --git a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp index dbfb9b1204d..ad1709bd933 100644 --- a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp +++ b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp @@ -32,7 +32,8 @@ TEST(MathMixMatFun, logSumExp) { stan::test::expect_ad(tols, f, x); Eigen::RowVectorXd rx = x; stan::test::expect_ad(tols, f, rx); - std::vector stx = std::vector(x.data(), x.data() + x.size()); + std::vector stx = std::vector(x.data(), + x.data() + x.size()); stan::test::expect_ad(tols, f, stx); } @@ -47,8 +48,10 @@ TEST(MathMixMatFun, logSumExp) { x2b.transpose(), x2c.transpose()}; stan::test::expect_ad(tols, f, strx); - std::vector> ststx{std::vector(x2.data(), x2.data() + x2.size()), - std::vector(x2b.data(), x2b.data() + x2b.size()), - std::vector(x2c.data(), x2c.data() + x2c.size())}; + std::vector> ststx{ + std::vector(x2.data(), x2.data() + x2.size()), + std::vector(x2b.data(), x2b.data() + x2b.size()), + std::vector(x2c.data(), x2c.data() + x2c.size()) + }; stan::test::expect_ad(tols, f, ststx); } diff --git a/test/unit/math/prim/mat/fun/head_test.cpp b/test/unit/math/prim/mat/fun/head_test.cpp index 53b56baa817..f95197bc6cb 100644 --- a/test/unit/math/prim/mat/fun/head_test.cpp +++ b/test/unit/math/prim/mat/fun/head_test.cpp @@ -54,7 +54,7 @@ TEST(MathMatrixHead, HeadRowVector4) { using stan::math::head; Eigen::RowVectorXd v(3); v << 1, 2, 3; - std::vector vind{1,2,1}; + std::vector vind{1, 2, 1}; std::vector st_v{v, v, v}; std::vector st_t = head(st_v, vind); diff --git a/test/unit/math/prim/mat/fun/log_softmax_test.cpp b/test/unit/math/prim/mat/fun/log_softmax_test.cpp index 2d28873958b..70b38902394 100644 --- a/test/unit/math/prim/mat/fun/log_softmax_test.cpp +++ b/test/unit/math/prim/mat/fun/log_softmax_test.cpp @@ -1,5 +1,6 @@ #include #include +#include void test_log_softmax(const Eigen::Matrix& theta) { using Eigen::Dynamic; @@ -33,32 +34,32 @@ TEST(MathMatrixPrimMat, log_softmax) { // x << 0.0; // test_log_softmax(x); - std::vector in{-1,1}; + std::vector in{-1, 1}; std::vector out = log_softmax(in); stan::math::vector_d x2(2); x2 << -1.0, 1.0; - stan::math::matrix_d m2(2,2); + stan::math::matrix_d m2(2, 2); m2 << -1.0, 1.0, -1.0, 1.0; stan::math::vector_d x2_out = log_softmax(x2); - EXPECT_FLOAT_EQ(out[0],x2_out[0]); - EXPECT_FLOAT_EQ(out[1],x2_out[1]); + EXPECT_FLOAT_EQ(out[0], x2_out[0]); + EXPECT_FLOAT_EQ(out[1], x2_out[1]); x2_out = log_softmax(m2.diagonal()); - EXPECT_FLOAT_EQ(out[0],x2_out[0]); - EXPECT_FLOAT_EQ(out[1],x2_out[1]); + EXPECT_FLOAT_EQ(out[0], x2_out[0]); + EXPECT_FLOAT_EQ(out[1], x2_out[1]); std::vector invec{x2, x2}; std::vector outvec = log_softmax(invec); std::vector> instvec{in, in}; std::vector> outstvec = log_softmax(instvec); - EXPECT_FLOAT_EQ(outvec[0][0],outstvec[0][0]); - EXPECT_FLOAT_EQ(outvec[0][1],outstvec[0][1]); - EXPECT_FLOAT_EQ(outvec[1][0],outstvec[1][0]); - EXPECT_FLOAT_EQ(outvec[1][1],outstvec[1][1]); + EXPECT_FLOAT_EQ(outvec[0][0], outstvec[0][0]); + EXPECT_FLOAT_EQ(outvec[0][1], outstvec[0][1]); + EXPECT_FLOAT_EQ(outvec[1][0], outstvec[1][0]); + EXPECT_FLOAT_EQ(outvec[1][1], outstvec[1][1]); // stan::math::vector_d x3(3); // x3 << -1.0, 1.0, 10.0; diff --git a/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp b/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp index 4ee758f5f0d..99eaba0bc55 100644 --- a/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp +++ b/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp @@ -1,6 +1,7 @@ #include #include #include +#include template void test_log_sum_exp(const Eigen::Matrix& as) { @@ -30,10 +31,9 @@ TEST(MathFunctions, log_sum_exp_mat) { double st2_out = log_sum_exp(st2); EXPECT_FLOAT_EQ(m1_out, st1_out); EXPECT_FLOAT_EQ(m2_out, st2_out); - //test_log_sum_exp(m); - std::vector st_m{m1,m2}; - std::vector> st_st{st1,st2}; + std::vector st_m{m1, m2}; + std::vector> st_st{st1, st2}; std::vector m_out = log_sum_exp(st_m); std::vector st_out = log_sum_exp(st_st); EXPECT_FLOAT_EQ(m1_out, m_out[0]); From d97c9633c7f7b573c5ff26b786d248b440d3a0ab Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Sun, 29 Dec 2019 10:20:16 +0800 Subject: [PATCH 08/22] Tidy missing doc --- stan/math/prim/mat.hpp | 1 + stan/math/prim/mat/err/check_column_index.hpp | 4 +- stan/math/prim/mat/err/check_row_index.hpp | 4 +- stan/math/prim/mat/fun/log_softmax.hpp | 2 +- .../prim/mat/vectorize/apply_vector_unary.hpp | 178 +++++++++++++----- stan/math/rev/core/matrix_vari.hpp | 7 +- test/unit/math/prim/mat/fun/head_test.cpp | 4 + .../math/prim/scal/fun/log_sum_exp_test.cpp | 2 +- 8 files changed, 144 insertions(+), 58 deletions(-) diff --git a/stan/math/prim/mat.hpp b/stan/math/prim/mat.hpp index 7acc453ef7e..51bffe100f0 100644 --- a/stan/math/prim/mat.hpp +++ b/stan/math/prim/mat.hpp @@ -333,6 +333,7 @@ #include #include +#include #include diff --git a/stan/math/prim/mat/err/check_column_index.hpp b/stan/math/prim/mat/err/check_column_index.hpp index 41d04212b97..4c6a0b29ae3 100644 --- a/stan/math/prim/mat/err/check_column_index.hpp +++ b/stan/math/prim/mat/err/check_column_index.hpp @@ -17,9 +17,7 @@ namespace math { * stan::error_index::value. This function will * throw an std::out_of_range exception if * the index is out of bounds. - * @tparam T_y Type of scalar - * @tparam R number of rows or Eigen::Dynamic - * @tparam C number of columns or Eigen::Dynamic + * @tparam Derived Type of Eigen object * @param function Function name (for error messages) * @param name Variable name (for error messages) * @param y matrix to test diff --git a/stan/math/prim/mat/err/check_row_index.hpp b/stan/math/prim/mat/err/check_row_index.hpp index c0d2ec006f1..616a4040f89 100644 --- a/stan/math/prim/mat/err/check_row_index.hpp +++ b/stan/math/prim/mat/err/check_row_index.hpp @@ -14,9 +14,7 @@ namespace math { * Check if the specified index is a valid row of the matrix * This check is 1-indexed by default. This behavior can be changed * by setting stan::error_index::value. - * @tparam T Scalar type - * @tparam R number of rows or Eigen::Dynamic - * @tparam C number of columns or Eigen::Dynamic + * @tparam Derived Type of Eigen object * @param function Function name (for error messages) * @param name Variable name (for error messages) * @param y matrix to test diff --git a/stan/math/prim/mat/fun/log_softmax.hpp b/stan/math/prim/mat/fun/log_softmax.hpp index 91195d3ddd5..9d4bc1d585c 100644 --- a/stan/math/prim/mat/fun/log_softmax.hpp +++ b/stan/math/prim/mat/fun/log_softmax.hpp @@ -35,7 +35,7 @@ namespace math { * \f$ * * @tparam T Type of input vector to transform. - * @param[in] v Vector to transform. + * @param[in] x Vector to transform. * @return Unit simplex result of the softmax transform of the vector. */ template >>...> diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp index ee036196df0..f3ae81de7c2 100644 --- a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -15,7 +15,8 @@ struct apply_vector_unary {}; /** * Base template class for vectorization of unary vector functions * defined by applying a functor to a standard library vector, Eigen dense - * matrix expression template, or container of these. + * matrix expression template, or container of these. For each specialisation, + * the same vector type as the input is returned. * * Three taxonomies of unary vector functions are implemented: * - f(vector) -> vector @@ -81,11 +82,12 @@ struct apply_vector_unary> { /** * Specialisation for use with (non-nested) std::vectors. Inputs are mapped * to Eigen column vectors and then passed to the base (Eigen) template. + * An std::vector (or scalar) is then returned as the result. */ template struct apply_vector_unary> { - using scalar_t = scalar_type_t; - using map_t = typename Eigen::Map>; + using T_scalar = scalar_type_t; + using T_map = typename Eigen::Map>; /** * Member function for applying a functor to a vector and subsequently @@ -98,17 +100,16 @@ struct apply_vector_unary> { * @return std::vector with result of applying functor to input. */ template - static inline std::vector apply(const T& x, const F& f) { - Eigen::Matrix result - = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); - return std::vector(result.data(), + static inline std::vector apply(const T& x, const F& f) { + Eigen::Matrix result + = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); + return std::vector(result.data(), result.data() + result.size()); } /** * Member function for applying a functor to a vector and a scalar - * and subsequently returning a vector. The 'auto' return type is - * used here so that an expression template is returned. + * and subsequently returning a vector. * * @tparam T Type of argument to which functor is applied. * @tparam T2 Type of scalar to pass to functor. @@ -119,13 +120,13 @@ struct apply_vector_unary> { * @return std::vector with result of applying functor and scalar to input. */ template - static inline std::vector apply_scalar(const T& x, + static inline std::vector apply_scalar(const T& x, const T2& y, const F& f) { - Eigen::Matrix result - = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), + Eigen::Matrix result + = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), y, f); - return std::vector(result.data(), + return std::vector(result.data(), result.data() + result.size()); } @@ -140,93 +141,176 @@ struct apply_vector_unary> { * @return scalar result of applying functor to input vector. */ template - static inline scalar_t reduce(const T& x, const F& f) { - return apply_vector_unary::reduce(as_column_vector_or_scalar(x), f); + static inline T_scalar reduce(const T& x, const F& f) { + return apply_vector_unary::reduce(as_column_vector_or_scalar(x), f); } }; /** - * Specialisation for use with std::vectors of Eigen types + * Specialisation for use with nested containers (std::vectors) of Eigen types. + * For each of the member functions, an std::vector with the appropriate + * type (vector or scalar) is returned. + * */ template struct apply_vector_unary> { - using eigen_t = value_type_t; - using scalar_t = typename eigen_t::Scalar; - using return_t = std::vector>; + using T_eigen = value_type_t; + using T_scalar = typename T_eigen::Scalar; + using T_return = std::vector>; + /** + * Member function for applying a functor to each Eigen type in an std::vector + * and subsequently returning an std::vector of evaluated Eigen expressions. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x std::vector of Eigen inputs to which operation is applied. + * @param f functor to apply to vector input. + * @return std::vector of Eigen objects with result of applying functor to + * input. + */ template - static inline return_t apply(const T& x, const F& f) { + static inline T_return apply(const T& x, const F& f) { size_t x_size = x.size(); - return_t result(x_size); + T_return result(x_size); for (size_t i = 0; i < x_size; ++i) - result[i] = apply_vector_unary::apply(x[i], f); + result[i] = apply_vector_unary::apply(x[i], f); return result; } + /** + * Member function for applying a functor to each Eigen type in an + * std::vector, as well as a scalar and subsequently returning an + * std::vector of evaluated Eigen expressions. The provided scalar can either + * be a single scalar which is applied to every Eigen object, or a vector of + * scalars to be applied, one for each Eigen object in the std::vector + * + * @tparam T Type of argument to which functor is applied. + * @tparam T2 Type of scalar (or vector) to pass to functor. + * @tparam F Type of functor to apply. + * @param x std::vector of Eigen inputs to which operation is applied. + * @param y scalar (or vector) passed to functor. + * @param f functor to apply to vector input. + * @return std::vector of Eigen objects with result of applying functor and + * scalar to input. + */ template - static inline return_t apply_scalar(const T& x, const T2& y, const F& f) { + static inline T_return apply_scalar(const T& x, const T2& y, const F& f) { scalar_seq_view y_vec(y); size_t x_size = x.size(); - return_t result(x_size); + T_return result(x_size); for (size_t i = 0; i < x_size; ++i) - result[i] = apply_vector_unary::apply_scalar(x[i], y_vec[i], f); + result[i] = apply_vector_unary::apply_scalar(x[i], y_vec[i], f); return result; } + /** + * Member function for applying a functor to each Eigen type in an + * std::vector and subsequently returning an std::vector of scalars. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x std::vector of Eigen inputs to which operation is applied. + * @param f functor to apply to vector input. + * @return std::vector of scalars with result of applying functor to input. + */ template - static inline std::vector reduce(const T& x, const F& f) { + static inline std::vector reduce(const T& x, const F& f) { size_t x_size = x.size(); - std::vector result(x_size); + std::vector result(x_size); for (size_t i = 0; i < x_size; ++i) - result[i] = apply_vector_unary::reduce(x[i], f); + result[i] = apply_vector_unary::reduce(x[i], f); return result; } }; +/** + * Specialisation for use with nested containers (std::vectors) of std::vectors. + * For each of the member functions, an std::vector with the appropriate + * type (vector or scalar) is returned. + * + */ template struct apply_vector_unary> { - using scalar_t = scalar_type_t; - using return_t = typename std::vector>; - using map_t = - typename Eigen::Map>; + using T_scalar = scalar_type_t; + using T_return = typename std::vector>; + using T_map = + typename Eigen::Map>; + /** + * Member function for applying a functor to each std::vector in an + * std::vector and subsequently returning an std::vector of std::vectors. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x std::vector of std::vector to which operation is applied. + * @param f functor to apply to vector input. + * @return std::vector of std::vectors with result of applying functor to + * input. + */ template - static inline return_t apply(const T& x, const F& f) { + static inline T_return apply(const T& x, const F& f) { size_t x_size = x.size(); - return_t result(x_size); - Eigen::Matrix inter; + T_return result(x_size); + Eigen::Matrix inter; for (size_t i = 0; i < x_size; ++i) { - inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), + inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), f); - result[i] = std::vector(inter.data(), + result[i] = std::vector(inter.data(), inter.data() + inter.size()); } return result; } + /** + * Member function for applying a functor to each std::vector in an + * std::vector, as well as a scalar and subsequently returning an + * std::vector of evaluated std::vectors. The provided scalar can either + * be a single scalar which is applied to every std::vector, or a vector of + * scalars to be applied, one for each std::vector in the std::vector + * + * @tparam T Type of argument to which functor is applied. + * @tparam T2 Type of scalar (or vector) to pass to functor. + * @tparam F Type of functor to apply. + * @param x std::vector of std::vectors to which operation is applied. + * @param y scalar (or vector) passed to functor. + * @param f functor to apply to vector input. + * @return std::vector of std::vectors with result of applying functor and + * scalar to input. + */ template - static inline return_t apply_scalar(const T& x, const T2& y, const F& f) { + static inline T_return apply_scalar(const T& x, const T2& y, const F& f) { scalar_seq_view y_vec(y); size_t x_size = x.size(); - return_t result(x_size); - Eigen::Matrix inter; + T_return result(x_size); + Eigen::Matrix inter; for (size_t i = 0; i < x_size; ++i) { - inter = apply_vector_unary::apply_scalar( + inter = apply_vector_unary::apply_scalar( as_column_vector_or_scalar(x[i]), y_vec[i], f); - result[i] = std::vector(inter.data(), + result[i] = std::vector(inter.data(), inter.data() + inter.size()); } return result; } + /** + * Member function for applying a functor to each std::vector in an + * std::vector and subsequently returning an std::vector of scalars. + * + * @tparam T Type of argument to which functor is applied. + * @tparam F Type of functor to apply. + * @param x std::vector of std::vectors to which operation is applied. + * @param f functor to apply to vector input. + * @return std::vector of scalars with result of applying functor to input. + */ template - static inline std::vector reduce(const T& x, const F& f) { + static inline std::vector reduce(const T& x, const F& f) { size_t x_size = x.size(); - std::vector result(x_size); + std::vector result(x_size); for (size_t i = 0; i < x_size; ++i) { - result[i] = apply_vector_unary::reduce( + result[i] = apply_vector_unary::reduce( as_column_vector_or_scalar(x[i]), f); } return result; diff --git a/stan/math/rev/core/matrix_vari.hpp b/stan/math/rev/core/matrix_vari.hpp index 95f41e99c7d..ad162dc8a2f 100644 --- a/stan/math/rev/core/matrix_vari.hpp +++ b/stan/math/rev/core/matrix_vari.hpp @@ -17,11 +17,12 @@ class op_matrix_vari : public vari { vari** vis_; public: - template >>...> - op_matrix_vari(double f, T&& vs) + template >...> + op_matrix_vari(double f, const Eigen::MatrixBase& vs) : vari(f), size_(vs.size()) { vis_ = ChainableStack::instance_->memalloc_.alloc_array(size_); - Eigen::Map(vis_, vs.rows(), vs.cols()) = vs.vi(); + Eigen::Map>(vis_, vs.rows(), + vs.cols()) = vs.vi(); } vari* operator[](size_t n) const { return vis_[n]; } size_t size() { return size_; } diff --git a/test/unit/math/prim/mat/fun/head_test.cpp b/test/unit/math/prim/mat/fun/head_test.cpp index f95197bc6cb..86b8fdd78e5 100644 --- a/test/unit/math/prim/mat/fun/head_test.cpp +++ b/test/unit/math/prim/mat/fun/head_test.cpp @@ -55,9 +55,11 @@ TEST(MathMatrixHead, HeadRowVector4) { Eigen::RowVectorXd v(3); v << 1, 2, 3; std::vector vind{1, 2, 1}; + std::vector ns{2, 2, 2}; std::vector st_v{v, v, v}; std::vector st_t = head(st_v, vind); + std::vector st_t2 = head(st_v, ns); Eigen::RowVectorXd v01 = head(v, 2); EXPECT_EQ(2, v01.size()); @@ -87,8 +89,10 @@ TEST(MathMatrixHead, HeadStdVector3) { v.push_back(1); v.push_back(2); v.push_back(3); + std::vector ns{2, 2, 2}; std::vector> st_v{v, v, v}; std::vector> st_t = head(st_v, 2); + std::vector> st_t2 = head(st_v, ns); EXPECT_THROW(head(v, 4), std::out_of_range); } TEST(MathMatrixHead, HeadStdVector4) { diff --git a/test/unit/math/prim/scal/fun/log_sum_exp_test.cpp b/test/unit/math/prim/scal/fun/log_sum_exp_test.cpp index 439c77e1a9b..9a59ce30de9 100644 --- a/test/unit/math/prim/scal/fun/log_sum_exp_test.cpp +++ b/test/unit/math/prim/scal/fun/log_sum_exp_test.cpp @@ -1,5 +1,5 @@ #include -#include +#include #include #include #include From 461f0b0f88fb3b2bb995588ff0d1c57f949398bc Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Sun, 29 Dec 2019 15:27:52 +0800 Subject: [PATCH 09/22] log_softmax doc errors --- stan/math/fwd/mat/fun/log_softmax.hpp | 3 +-- stan/math/prim/mat/fun/log_softmax.hpp | 2 +- stan/math/rev/mat/fun/log_softmax.hpp | 3 +-- 3 files changed, 3 insertions(+), 5 deletions(-) diff --git a/stan/math/fwd/mat/fun/log_softmax.hpp b/stan/math/fwd/mat/fun/log_softmax.hpp index b40887776bf..79a3284c6cf 100644 --- a/stan/math/fwd/mat/fun/log_softmax.hpp +++ b/stan/math/fwd/mat/fun/log_softmax.hpp @@ -12,8 +12,7 @@ namespace stan { namespace math { /** - * Return the softmax of the specified vector or container of vectors. - * Softmax is guaranteed to return a simplex. + * Return the log softmax of the specified vector or container of vectors. * * @tparam T Type of input vector or matrix. * @param[in] x Unconstrained input vector. diff --git a/stan/math/prim/mat/fun/log_softmax.hpp b/stan/math/prim/mat/fun/log_softmax.hpp index 9d4bc1d585c..ad2829551b5 100644 --- a/stan/math/prim/mat/fun/log_softmax.hpp +++ b/stan/math/prim/mat/fun/log_softmax.hpp @@ -36,7 +36,7 @@ namespace math { * * @tparam T Type of input vector to transform. * @param[in] x Vector to transform. - * @return Unit simplex result of the softmax transform of the vector. + * @return log unit simplex result of the softmax transform of the vector. */ template >>...> inline auto log_softmax(T&& x) { diff --git a/stan/math/rev/mat/fun/log_softmax.hpp b/stan/math/rev/mat/fun/log_softmax.hpp index 64834de07ff..370fa54bd3a 100644 --- a/stan/math/rev/mat/fun/log_softmax.hpp +++ b/stan/math/rev/mat/fun/log_softmax.hpp @@ -46,8 +46,7 @@ class log_softmax_elt_vari : public vari { } // namespace internal /** - * Return the softmax of the specified Eigen vector. Softmax is - * guaranteed to return a simplex. + * Return the log softmax of the specified vector or container of vectors. * * The gradient calculations are unfolded. * From 9fa51242de44e01b2ad7f409e1ae3e222f264fa0 Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Sun, 29 Dec 2019 07:51:55 +0000 Subject: [PATCH 10/22] [Jenkins] auto-formatting by clang-format version 5.0.0-3~16.04.1 (tags/RELEASE_500/final) --- stan/math/fwd/mat/fun/log_softmax.hpp | 2 +- stan/math/fwd/mat/fun/log_sum_exp.hpp | 6 +-- stan/math/prim/mat/fun/head.hpp | 22 +++++----- stan/math/prim/mat/fun/log_softmax.hpp | 3 +- stan/math/prim/mat/fun/log_sum_exp.hpp | 2 +- .../prim/mat/vectorize/apply_vector_unary.hpp | 44 +++++++++---------- stan/math/rev/core/matrix_vari.hpp | 4 +- stan/math/rev/mat/fun/log_softmax.hpp | 4 +- stan/math/rev/mat/fun/log_sum_exp.hpp | 2 +- .../math/mix/mat/fun/log_sum_exp_test.cpp | 15 +++---- .../math/prim/mat/fun/log_sum_exp_test.cpp | 2 +- 11 files changed, 49 insertions(+), 57 deletions(-) diff --git a/stan/math/fwd/mat/fun/log_softmax.hpp b/stan/math/fwd/mat/fun/log_softmax.hpp index 79a3284c6cf..a2f474fe86d 100644 --- a/stan/math/fwd/mat/fun/log_softmax.hpp +++ b/stan/math/fwd/mat/fun/log_softmax.hpp @@ -21,7 +21,7 @@ namespace math { */ template >>...> inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [&](auto& alpha){ + return apply_vector_unary::apply(std::forward(x), [&](auto& alpha) { using T_fvar = value_type_t; using T_fvar_inner = typename T_fvar::Scalar; diff --git a/stan/math/fwd/mat/fun/log_sum_exp.hpp b/stan/math/fwd/mat/fun/log_sum_exp.hpp index 6eb936a9804..3245dc10113 100644 --- a/stan/math/fwd/mat/fun/log_sum_exp.hpp +++ b/stan/math/fwd/mat/fun/log_sum_exp.hpp @@ -26,14 +26,14 @@ namespace math { */ template >>...> inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [&](auto& v){ + return apply_vector_unary::reduce(std::forward(x), [&](auto& v) { using T_fvar_inner = typename value_type_t::Scalar; using mat_type = Eigen::Matrix; mat_type vals = v.val(); mat_type exp_vals = vals.array().exp(); - return fvar(log_sum_exp(vals), - v.d().cwiseProduct(exp_vals).sum() / exp_vals.sum()); + return fvar( + log_sum_exp(vals), v.d().cwiseProduct(exp_vals).sum() / exp_vals.sum()); }); } diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index d32f429a5ec..c8b429c7cf9 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -24,17 +24,17 @@ namespace math { */ template inline auto head(T&& x, const T2& n) { - return apply_vector_unary::apply_scalar(std::forward(x), n, - [&](auto& v, auto& m){ - if (m != 0) { - if (v.rows() == 1) { - check_column_index("head", "n", v, m); - } else { - check_row_index("head", "n", v, m); - } - } - return v.head(m); - }); + return apply_vector_unary::apply_scalar( + std::forward(x), n, [&](auto& v, auto& m) { + if (m != 0) { + if (v.rows() == 1) { + check_column_index("head", "n", v, m); + } else { + check_row_index("head", "n", v, m); + } + } + return v.head(m); + }); } } // namespace math diff --git a/stan/math/prim/mat/fun/log_softmax.hpp b/stan/math/prim/mat/fun/log_softmax.hpp index ad2829551b5..0e111145594 100644 --- a/stan/math/prim/mat/fun/log_softmax.hpp +++ b/stan/math/prim/mat/fun/log_softmax.hpp @@ -6,7 +6,6 @@ #include #include - namespace stan { namespace math { @@ -40,7 +39,7 @@ namespace math { */ template >>...> inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [&](auto& v){ + return apply_vector_unary::apply(std::forward(x), [&](auto& v) { check_nonzero_size("log_softmax", "v", v); return (v.array() - log_sum_exp(v)).matrix(); }); diff --git a/stan/math/prim/mat/fun/log_sum_exp.hpp b/stan/math/prim/mat/fun/log_sum_exp.hpp index fe6d25ce013..02ac35d55d4 100644 --- a/stan/math/prim/mat/fun/log_sum_exp.hpp +++ b/stan/math/prim/mat/fun/log_sum_exp.hpp @@ -27,7 +27,7 @@ namespace math { */ template >>...> inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [&](auto& v){ + return apply_vector_unary::reduce(std::forward(x), [&](auto& v) { if (v.size() == 0) { return -std::numeric_limits::infinity(); } diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp index f3ae81de7c2..fed3d130627 100644 --- a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -91,7 +91,7 @@ struct apply_vector_unary> { /** * Member function for applying a functor to a vector and subsequently - * returning a vector. + * returning a vector. * * @tparam T Type of argument to which functor is applied. * @tparam F Type of functor to apply. @@ -102,14 +102,13 @@ struct apply_vector_unary> { template static inline std::vector apply(const T& x, const F& f) { Eigen::Matrix result - = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); - return std::vector(result.data(), - result.data() + result.size()); + = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); + return std::vector(result.data(), result.data() + result.size()); } /** * Member function for applying a functor to a vector and a scalar - * and subsequently returning a vector. + * and subsequently returning a vector. * * @tparam T Type of argument to which functor is applied. * @tparam T2 Type of scalar to pass to functor. @@ -120,14 +119,12 @@ struct apply_vector_unary> { * @return std::vector with result of applying functor and scalar to input. */ template - static inline std::vector apply_scalar(const T& x, - const T2& y, + static inline std::vector apply_scalar(const T& x, const T2& y, const F& f) { Eigen::Matrix result - = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), - y, f); - return std::vector(result.data(), - result.data() + result.size()); + = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), + y, f); + return std::vector(result.data(), result.data() + result.size()); } /** @@ -150,15 +147,15 @@ struct apply_vector_unary> { * Specialisation for use with nested containers (std::vectors) of Eigen types. * For each of the member functions, an std::vector with the appropriate * type (vector or scalar) is returned. - * + * */ template struct apply_vector_unary> { using T_eigen = value_type_t; using T_scalar = typename T_eigen::Scalar; - using T_return = std::vector>; + using T_return + = std::vector>; /** * Member function for applying a functor to each Eigen type in an std::vector @@ -230,14 +227,13 @@ struct apply_vector_unary> { * Specialisation for use with nested containers (std::vectors) of std::vectors. * For each of the member functions, an std::vector with the appropriate * type (vector or scalar) is returned. - * + * */ template struct apply_vector_unary> { using T_scalar = scalar_type_t; using T_return = typename std::vector>; - using T_map = - typename Eigen::Map>; + using T_map = typename Eigen::Map>; /** * Member function for applying a functor to each std::vector in an @@ -258,8 +254,8 @@ struct apply_vector_unary> { for (size_t i = 0; i < x_size; ++i) { inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), f); - result[i] = std::vector(inter.data(), - inter.data() + inter.size()); + result[i] + = std::vector(inter.data(), inter.data() + inter.size()); } return result; } @@ -288,9 +284,9 @@ struct apply_vector_unary> { Eigen::Matrix inter; for (size_t i = 0; i < x_size; ++i) { inter = apply_vector_unary::apply_scalar( - as_column_vector_or_scalar(x[i]), y_vec[i], f); - result[i] = std::vector(inter.data(), - inter.data() + inter.size()); + as_column_vector_or_scalar(x[i]), y_vec[i], f); + result[i] + = std::vector(inter.data(), inter.data() + inter.size()); } return result; } @@ -311,7 +307,7 @@ struct apply_vector_unary> { std::vector result(x_size); for (size_t i = 0; i < x_size; ++i) { result[i] = apply_vector_unary::reduce( - as_column_vector_or_scalar(x[i]), f); + as_column_vector_or_scalar(x[i]), f); } return result; } diff --git a/stan/math/rev/core/matrix_vari.hpp b/stan/math/rev/core/matrix_vari.hpp index ad162dc8a2f..3b47b721599 100644 --- a/stan/math/rev/core/matrix_vari.hpp +++ b/stan/math/rev/core/matrix_vari.hpp @@ -21,8 +21,8 @@ class op_matrix_vari : public vari { op_matrix_vari(double f, const Eigen::MatrixBase& vs) : vari(f), size_(vs.size()) { vis_ = ChainableStack::instance_->memalloc_.alloc_array(size_); - Eigen::Map>(vis_, vs.rows(), - vs.cols()) = vs.vi(); + Eigen::Map>(vis_, vs.rows(), vs.cols()) + = vs.vi(); } vari* operator[](size_t n) const { return vis_[n]; } size_t size() { return size_; } diff --git a/stan/math/rev/mat/fun/log_softmax.hpp b/stan/math/rev/mat/fun/log_softmax.hpp index 370fa54bd3a..d19969dd4b9 100644 --- a/stan/math/rev/mat/fun/log_softmax.hpp +++ b/stan/math/rev/mat/fun/log_softmax.hpp @@ -57,7 +57,7 @@ class log_softmax_elt_vari : public vari { */ template >>...> inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [&](auto& alpha){ + return apply_vector_unary::apply(std::forward(x), [&](auto& alpha) { const int a_size = alpha.size(); check_nonzero_size("log_softmax", "alpha", alpha); @@ -80,7 +80,7 @@ inline auto log_softmax(T&& x) { double* softmax_alpha_d_array = ChainableStack::instance_->memalloc_.alloc_array(a_size); Eigen::Map(softmax_alpha_d_array, a_size) - = softmax_alpha_d.array() / sum; + = softmax_alpha_d.array() / sum; vector_v log_softmax_alpha(a_size); for (int k = 0; k < a_size; ++k) { diff --git a/stan/math/rev/mat/fun/log_sum_exp.hpp b/stan/math/rev/mat/fun/log_sum_exp.hpp index fcaf8da39af..686afca8dc1 100644 --- a/stan/math/rev/mat/fun/log_sum_exp.hpp +++ b/stan/math/rev/mat/fun/log_sum_exp.hpp @@ -36,7 +36,7 @@ class log_sum_exp_matrix_vari : public op_matrix_vari { */ template >>...> inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [&](auto& v){ + return apply_vector_unary::reduce(std::forward(x), [&](auto& v) { return var(new internal::log_sum_exp_matrix_vari(v)); }); } diff --git a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp index ad1709bd933..c22d20a3aed 100644 --- a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp +++ b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp @@ -32,26 +32,23 @@ TEST(MathMixMatFun, logSumExp) { stan::test::expect_ad(tols, f, x); Eigen::RowVectorXd rx = x; stan::test::expect_ad(tols, f, rx); - std::vector stx = std::vector(x.data(), - x.data() + x.size()); + std::vector stx + = std::vector(x.data(), x.data() + x.size()); stan::test::expect_ad(tols, f, stx); } - Eigen::MatrixXd x23(2, 2); x23 << 1, 2, 3, 4; stan::test::expect_ad(f, x23); std::vector stvx{x2, x2b, x2c}; stan::test::expect_ad(tols, f, stvx); - std::vector strx{x2.transpose(), - x2b.transpose(), + std::vector strx{x2.transpose(), x2b.transpose(), x2c.transpose()}; stan::test::expect_ad(tols, f, strx); std::vector> ststx{ - std::vector(x2.data(), x2.data() + x2.size()), - std::vector(x2b.data(), x2b.data() + x2b.size()), - std::vector(x2c.data(), x2c.data() + x2c.size()) - }; + std::vector(x2.data(), x2.data() + x2.size()), + std::vector(x2b.data(), x2b.data() + x2b.size()), + std::vector(x2c.data(), x2c.data() + x2c.size())}; stan::test::expect_ad(tols, f, ststx); } diff --git a/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp b/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp index 99eaba0bc55..21b14716498 100644 --- a/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp +++ b/test/unit/math/prim/mat/fun/log_sum_exp_test.cpp @@ -15,9 +15,9 @@ void test_log_sum_exp(const Eigen::Matrix& as) { } TEST(MathFunctions, log_sum_exp_mat) { - using stan::math::log_sum_exp; using Eigen::Dynamic; using Eigen::Matrix; + using stan::math::log_sum_exp; Matrix m1(3, 2); m1 << 1, 2, 3, 4, 5, 6; From cba3a3055157c439b5d5be441fc3c8307efc44c6 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Mon, 30 Dec 2019 10:13:13 +0800 Subject: [PATCH 11/22] Fix failing test --- stan/math/prim/mat/fun/head.hpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index c8b429c7cf9..887866de33d 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -23,9 +23,9 @@ namespace math { * @throw std::out_of_range if n is out of range. */ template -inline auto head(T&& x, const T2& n) { +inline auto head(T&& x, T2&& n) { return apply_vector_unary::apply_scalar( - std::forward(x), n, [&](auto& v, auto& m) { + std::forward(x), std::forward(n), [&](auto& v, auto& m) { if (m != 0) { if (v.rows() == 1) { check_column_index("head", "n", v, m); @@ -33,7 +33,7 @@ inline auto head(T&& x, const T2& n) { check_row_index("head", "n", v, m); } } - return v.head(m); + return v.head(m).eval(); }); } From 5a934fe179a4a3f0fa2d889e907b11b41b599d66 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 31 Dec 2019 16:01:06 +0800 Subject: [PATCH 12/22] Revert head replacement --- stan/math/prim/mat/err/check_column_index.hpp | 8 +- stan/math/prim/mat/err/check_row_index.hpp | 8 +- stan/math/prim/mat/fun/head.hpp | 70 ++++-- .../prim/mat/vectorize/apply_vector_unary.hpp | 205 ++---------------- test/unit/math/mix/mat/fun/head_test.cpp | 12 +- test/unit/math/prim/mat/fun/head_test.cpp | 11 - 6 files changed, 87 insertions(+), 227 deletions(-) diff --git a/stan/math/prim/mat/err/check_column_index.hpp b/stan/math/prim/mat/err/check_column_index.hpp index 4c6a0b29ae3..13a8127765d 100644 --- a/stan/math/prim/mat/err/check_column_index.hpp +++ b/stan/math/prim/mat/err/check_column_index.hpp @@ -17,16 +17,18 @@ namespace math { * stan::error_index::value. This function will * throw an std::out_of_range exception if * the index is out of bounds. - * @tparam Derived Type of Eigen object + * @tparam T_y Type of scalar + * @tparam R number of rows or Eigen::Dynamic + * @tparam C number of columns or Eigen::Dynamic * @param function Function name (for error messages) * @param name Variable name (for error messages) * @param y matrix to test * @param i column index to check * @throw std::out_of_range if index is an invalid column */ -template +template inline void check_column_index(const char* function, const char* name, - const Eigen::MatrixBase& y, size_t i) { + const Eigen::Matrix& y, size_t i) { if (i >= stan::error_index::value && i < static_cast(y.cols()) + stan::error_index::value) { return; diff --git a/stan/math/prim/mat/err/check_row_index.hpp b/stan/math/prim/mat/err/check_row_index.hpp index 616a4040f89..7df742b7cc7 100644 --- a/stan/math/prim/mat/err/check_row_index.hpp +++ b/stan/math/prim/mat/err/check_row_index.hpp @@ -14,16 +14,18 @@ namespace math { * Check if the specified index is a valid row of the matrix * This check is 1-indexed by default. This behavior can be changed * by setting stan::error_index::value. - * @tparam Derived Type of Eigen object + * @tparam T Scalar type + * @tparam R number of rows or Eigen::Dynamic + * @tparam C number of columns or Eigen::Dynamic * @param function Function name (for error messages) * @param name Variable name (for error messages) * @param y matrix to test * @param i row index to check * @throw std::out_of_range if the index is out of range. */ -template +template inline void check_row_index(const char* function, const char* name, - const Eigen::MatrixBase& y, size_t i) { + const Eigen::Matrix& y, size_t i) { if (i >= stan::error_index::value && i < static_cast(y.rows()) + stan::error_index::value) { return; diff --git a/stan/math/prim/mat/fun/head.hpp b/stan/math/prim/mat/fun/head.hpp index 85219b23db4..fe60f072347 100644 --- a/stan/math/prim/mat/fun/head.hpp +++ b/stan/math/prim/mat/fun/head.hpp @@ -3,7 +3,7 @@ #include #include -#include +#include namespace stan { namespace math { @@ -12,28 +12,64 @@ namespace math { * Return the specified number of elements as a vector * from the front of the specified vector. * - * @tparam T Type of input vector. - * @tparam T2 Type of size variable. - * @param x Vector input. + * @tparam T type of elements in the vector + * @param v Vector input. * @param n Size of return. * @return The first n elements of v. * @throw std::out_of_range if n is out of range. */ -template -inline auto head(T&& x, T2&& n) { - return apply_vector_unary::apply_scalar( - std::forward(x), std::forward(n), [&](auto& v, auto& m) { - if (m != 0) { - if (v.rows() == 1) { - check_column_index("head", "n", v, m); - } else { - check_row_index("head", "n", v, m); - } - } - return v.head(m).eval(); - }); +template +inline Eigen::Matrix head( + const Eigen::Matrix& v, size_t n) { + if (n != 0) { + check_row_index("head", "n", v, n); + } + return v.head(n); +} + +/** + * Return the specified number of elements as a row vector + * from the front of the specified row vector. + * + * @tparam T type of elements in the vector + * @param rv Row vector. + * @param n Size of return row vector. + * @return The first n elements of rv. + * @throw std::out_of_range if n is out of range. + */ +template +inline Eigen::Matrix head( + const Eigen::Matrix& rv, size_t n) { + if (n != 0) { + check_column_index("head", "n", rv, n); + } + return rv.head(n); +} + +/** + * Return the specified number of elements as a standard vector + * from the front of the specified standard vector. + * + * @tparam T type of elements in the vector + * @param sv Standard vector. + * @param n Size of return. + * @return The first n elements of sv. + * @throw std::out_of_range if n is out of range. + */ +template +std::vector head(const std::vector& sv, size_t n) { + if (n != 0) { + check_std_vector_index("head", "n", sv, n); + } + + std::vector s; + for (size_t i = 0; i < n; ++i) { + s.push_back(sv[i]); + } + return s; } } // namespace math } // namespace stan + #endif diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp index fed3d130627..b7c0b55013b 100644 --- a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -18,10 +18,9 @@ struct apply_vector_unary {}; * matrix expression template, or container of these. For each specialisation, * the same vector type as the input is returned. * - * Three taxonomies of unary vector functions are implemented: + * Two taxonomies of unary vector functions are implemented: * - f(vector) -> vector * - f(vector) -> scalar - * - f(vector, scalar) -> vector * * This base template class takes (and returns) Eigen expression templates. */ @@ -44,25 +43,6 @@ struct apply_vector_unary> { return f(x); } - /** - * Member function for applying a functor to a vector and a scalar - * and subsequently returning a vector. The 'auto' return type is - * used here so that an expression template is returned. - * - * @tparam T Type of argument to which functor is applied. - * @tparam T2 Type of scalar to pass to functor. - * @tparam F Type of functor to apply. - * @param x Eigen input to which operation is applied. - * @param y scalar passed to functor. - * @param f functor to apply to Eigen input. - * @return Eigen expression template with result of applying functor - * and scalar to input - */ - template - static inline auto apply_scalar(const T& x, const T2& y, const F& f) { - return f(x, y); - } - /** * Member function for applying a functor to a vector and subsequently * returning a scalar. @@ -86,8 +66,8 @@ struct apply_vector_unary> { */ template struct apply_vector_unary> { - using T_scalar = scalar_type_t; - using T_map = typename Eigen::Map>; + using T_vt = value_type_t; + using T_map = typename Eigen::Map>; /** * Member function for applying a functor to a vector and subsequently @@ -100,31 +80,10 @@ struct apply_vector_unary> { * @return std::vector with result of applying functor to input. */ template - static inline std::vector apply(const T& x, const F& f) { - Eigen::Matrix result + static inline std::vector apply(const T& x, const F& f) { + Eigen::Matrix result = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); - return std::vector(result.data(), result.data() + result.size()); - } - - /** - * Member function for applying a functor to a vector and a scalar - * and subsequently returning a vector. - * - * @tparam T Type of argument to which functor is applied. - * @tparam T2 Type of scalar to pass to functor. - * @tparam F Type of functor to apply. - * @param x std::vector input to which operation is applied. - * @param y scalar passed to functor. - * @param f functor to apply to vector input. - * @return std::vector with result of applying functor and scalar to input. - */ - template - static inline std::vector apply_scalar(const T& x, const T2& y, - const F& f) { - Eigen::Matrix result - = apply_vector_unary::apply_scalar(as_column_vector_or_scalar(x), - y, f); - return std::vector(result.data(), result.data() + result.size()); + return std::vector(result.data(), result.data() + result.size()); } /** @@ -138,180 +97,62 @@ struct apply_vector_unary> { * @return scalar result of applying functor to input vector. */ template - static inline T_scalar reduce(const T& x, const F& f) { + static inline T_vt reduce(const T& x, const F& f) { return apply_vector_unary::reduce(as_column_vector_or_scalar(x), f); } }; /** - * Specialisation for use with nested containers (std::vectors) of Eigen types. + * Specialisation for use with nested containers (std::vectors). * For each of the member functions, an std::vector with the appropriate * type (vector or scalar) is returned. * */ template -struct apply_vector_unary> { - using T_eigen = value_type_t; - using T_scalar = typename T_eigen::Scalar; - using T_return - = std::vector>; +struct apply_vector_unary> { + using T_vt = value_type_t; + using T_st = value_type_t; /** - * Member function for applying a functor to each Eigen type in an std::vector - * and subsequently returning an std::vector of evaluated Eigen expressions. + * Member function for applying a functor to each container in an std::vector + * and subsequently returning an std::vector of containers. * * @tparam T Type of argument to which functor is applied. * @tparam F Type of functor to apply. - * @param x std::vector of Eigen inputs to which operation is applied. + * @param x std::vector of containers to which operation is applied. * @param f functor to apply to vector input. - * @return std::vector of Eigen objects with result of applying functor to + * @return std::vector of containers with result of applying functor to * input. */ template - static inline T_return apply(const T& x, const F& f) { - size_t x_size = x.size(); - T_return result(x_size); - for (size_t i = 0; i < x_size; ++i) - result[i] = apply_vector_unary::apply(x[i], f); - return result; - } - - /** - * Member function for applying a functor to each Eigen type in an - * std::vector, as well as a scalar and subsequently returning an - * std::vector of evaluated Eigen expressions. The provided scalar can either - * be a single scalar which is applied to every Eigen object, or a vector of - * scalars to be applied, one for each Eigen object in the std::vector - * - * @tparam T Type of argument to which functor is applied. - * @tparam T2 Type of scalar (or vector) to pass to functor. - * @tparam F Type of functor to apply. - * @param x std::vector of Eigen inputs to which operation is applied. - * @param y scalar (or vector) passed to functor. - * @param f functor to apply to vector input. - * @return std::vector of Eigen objects with result of applying functor and - * scalar to input. - */ - template - static inline T_return apply_scalar(const T& x, const T2& y, const F& f) { - scalar_seq_view y_vec(y); + static inline std::vector apply(const T& x, const F& f) { size_t x_size = x.size(); - T_return result(x_size); + std::vector result(x_size); for (size_t i = 0; i < x_size; ++i) - result[i] = apply_vector_unary::apply_scalar(x[i], y_vec[i], f); + result[i] = apply_vector_unary::apply(x[i], f); return result; } /** - * Member function for applying a functor to each Eigen type in an + * Member function for applying a functor to each container in an * std::vector and subsequently returning an std::vector of scalars. * * @tparam T Type of argument to which functor is applied. * @tparam F Type of functor to apply. - * @param x std::vector of Eigen inputs to which operation is applied. + * @param x std::vector of containers to which operation is applied. * @param f functor to apply to vector input. * @return std::vector of scalars with result of applying functor to input. */ template - static inline std::vector reduce(const T& x, const F& f) { + static inline std::vector reduce(const T& x, const F& f) { size_t x_size = x.size(); - std::vector result(x_size); + std::vector result(x_size); for (size_t i = 0; i < x_size; ++i) - result[i] = apply_vector_unary::reduce(x[i], f); + result[i] = apply_vector_unary::reduce(x[i], f); return result; } }; -/** - * Specialisation for use with nested containers (std::vectors) of std::vectors. - * For each of the member functions, an std::vector with the appropriate - * type (vector or scalar) is returned. - * - */ -template -struct apply_vector_unary> { - using T_scalar = scalar_type_t; - using T_return = typename std::vector>; - using T_map = typename Eigen::Map>; - - /** - * Member function for applying a functor to each std::vector in an - * std::vector and subsequently returning an std::vector of std::vectors. - * - * @tparam T Type of argument to which functor is applied. - * @tparam F Type of functor to apply. - * @param x std::vector of std::vector to which operation is applied. - * @param f functor to apply to vector input. - * @return std::vector of std::vectors with result of applying functor to - * input. - */ - template - static inline T_return apply(const T& x, const F& f) { - size_t x_size = x.size(); - T_return result(x_size); - Eigen::Matrix inter; - for (size_t i = 0; i < x_size; ++i) { - inter = apply_vector_unary::apply(as_column_vector_or_scalar(x[i]), - f); - result[i] - = std::vector(inter.data(), inter.data() + inter.size()); - } - return result; - } - - /** - * Member function for applying a functor to each std::vector in an - * std::vector, as well as a scalar and subsequently returning an - * std::vector of evaluated std::vectors. The provided scalar can either - * be a single scalar which is applied to every std::vector, or a vector of - * scalars to be applied, one for each std::vector in the std::vector - * - * @tparam T Type of argument to which functor is applied. - * @tparam T2 Type of scalar (or vector) to pass to functor. - * @tparam F Type of functor to apply. - * @param x std::vector of std::vectors to which operation is applied. - * @param y scalar (or vector) passed to functor. - * @param f functor to apply to vector input. - * @return std::vector of std::vectors with result of applying functor and - * scalar to input. - */ - template - static inline T_return apply_scalar(const T& x, const T2& y, const F& f) { - scalar_seq_view y_vec(y); - size_t x_size = x.size(); - T_return result(x_size); - Eigen::Matrix inter; - for (size_t i = 0; i < x_size; ++i) { - inter = apply_vector_unary::apply_scalar( - as_column_vector_or_scalar(x[i]), y_vec[i], f); - result[i] - = std::vector(inter.data(), inter.data() + inter.size()); - } - return result; - } - - /** - * Member function for applying a functor to each std::vector in an - * std::vector and subsequently returning an std::vector of scalars. - * - * @tparam T Type of argument to which functor is applied. - * @tparam F Type of functor to apply. - * @param x std::vector of std::vectors to which operation is applied. - * @param f functor to apply to vector input. - * @return std::vector of scalars with result of applying functor to input. - */ - template - static inline std::vector reduce(const T& x, const F& f) { - size_t x_size = x.size(); - std::vector result(x_size); - for (size_t i = 0; i < x_size; ++i) { - result[i] = apply_vector_unary::reduce( - as_column_vector_or_scalar(x[i]), f); - } - return result; - } -}; } // namespace math } // namespace stan diff --git a/test/unit/math/mix/mat/fun/head_test.cpp b/test/unit/math/mix/mat/fun/head_test.cpp index 893a880ad87..a5a4305d105 100644 --- a/test/unit/math/mix/mat/fun/head_test.cpp +++ b/test/unit/math/mix/mat/fun/head_test.cpp @@ -2,26 +2,16 @@ #include auto f(int i) { - return [=](const auto& y) { return stan::math::head(y, i).eval(); }; -} - -auto g(int i) { return [=](const auto& y) { return stan::math::head(y, i); }; } template void expect_head(const T& x, int n) { - std::vector> stx{x, x, x}; Eigen::VectorXd v = stan::test::to_vector(x); - std::vector stv{v, v, v}; Eigen::RowVectorXd rv = stan::test::to_row_vector(x); - std::vector strv{rv, rv, rv}; - stan::test::expect_ad(g(n), x); - stan::test::expect_ad(g(n), stx); + stan::test::expect_ad(f(n), x); stan::test::expect_ad(f(n), v); - stan::test::expect_ad(g(n), stv); stan::test::expect_ad(f(n), rv); - stan::test::expect_ad(g(n), strv); } TEST(MathMixMatFun, head) { diff --git a/test/unit/math/prim/mat/fun/head_test.cpp b/test/unit/math/prim/mat/fun/head_test.cpp index 86b8fdd78e5..6f5910e3291 100644 --- a/test/unit/math/prim/mat/fun/head_test.cpp +++ b/test/unit/math/prim/mat/fun/head_test.cpp @@ -8,7 +8,6 @@ TEST(MathMatrixHead, HeadVector1) { v << 1, 2, 3; EXPECT_EQ(0, head(v, 0).size()); } - TEST(MathMatrixHead, HeadVector2) { using stan::math::head; Eigen::VectorXd v(3); @@ -54,12 +53,6 @@ TEST(MathMatrixHead, HeadRowVector4) { using stan::math::head; Eigen::RowVectorXd v(3); v << 1, 2, 3; - std::vector vind{1, 2, 1}; - std::vector ns{2, 2, 2}; - - std::vector st_v{v, v, v}; - std::vector st_t = head(st_v, vind); - std::vector st_t2 = head(st_v, ns); Eigen::RowVectorXd v01 = head(v, 2); EXPECT_EQ(2, v01.size()); @@ -89,10 +82,6 @@ TEST(MathMatrixHead, HeadStdVector3) { v.push_back(1); v.push_back(2); v.push_back(3); - std::vector ns{2, 2, 2}; - std::vector> st_v{v, v, v}; - std::vector> st_t = head(st_v, 2); - std::vector> st_t2 = head(st_v, ns); EXPECT_THROW(head(v, 4), std::out_of_range); } TEST(MathMatrixHead, HeadStdVector4) { From 8ebea407dc6bea63097f2331cf0693a9da6f720e Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Tue, 31 Dec 2019 08:08:27 +0000 Subject: [PATCH 13/22] [Jenkins] auto-formatting by clang-format version 5.0.0-3~16.04.1 (tags/RELEASE_500/final) --- stan/math/prim/mat/vectorize/apply_vector_unary.hpp | 1 - 1 file changed, 1 deletion(-) diff --git a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp index b7c0b55013b..6ef78600196 100644 --- a/stan/math/prim/mat/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/mat/vectorize/apply_vector_unary.hpp @@ -153,7 +153,6 @@ struct apply_vector_unary> { } }; - } // namespace math } // namespace stan #endif From d776eac625e7cc240ce9648649f2ca3b7ea506f7 Mon Sep 17 00:00:00 2001 From: rok-cesnovar Date: Sun, 5 Jan 2020 11:35:49 +0100 Subject: [PATCH 14/22] merge conflicts fix --- stan/math/fwd/fun/log_softmax.hpp | 4 ++-- stan/math/fwd/fun/log_sum_exp.hpp | 31 +++++++++++++++++++++++++------ 2 files changed, 27 insertions(+), 8 deletions(-) diff --git a/stan/math/fwd/fun/log_softmax.hpp b/stan/math/fwd/fun/log_softmax.hpp index 2aa2fef9685..7ba2ec759fb 100644 --- a/stan/math/fwd/fun/log_softmax.hpp +++ b/stan/math/fwd/fun/log_softmax.hpp @@ -1,5 +1,5 @@ -#ifndef STAN_MATH_FWD_FUN_LOG_SOFTMAX_HPP -#define STAN_MATH_FWD_FUN_LOG_SOFTMAX_HPP +#ifndef STAN_MATH_FWD_LOG_SOFTMAX_HPP +#define STAN_MATH_FWD_LOG_SOFTMAX_HPP #include #include diff --git a/stan/math/fwd/fun/log_sum_exp.hpp b/stan/math/fwd/fun/log_sum_exp.hpp index 4c8ad3df687..742045e594a 100644 --- a/stan/math/fwd/fun/log_sum_exp.hpp +++ b/stan/math/fwd/fun/log_sum_exp.hpp @@ -57,13 +57,32 @@ fvar log_sum_exp(const std::vector >& v) { return fvar(log_sum_exp(vals), deriv / denominator); } -template -fvar log_sum_exp(const Eigen::Matrix, R, C>& v) { - Eigen::Matrix vals = v.val(); - Eigen::Matrix exp_vals = vals.array().exp(); +/** + * Return the log of the sum of the exponentiated values of the specified + * matrix of values. The matrix may be a full matrix, a vector, + * a row vector, or a container of these. + * + * The function is defined as follows to prevent overflow in exponential + * calculations. + * + * \f$\log \sum_{n=1}^N \exp(x_n) = \max(x) + \log \sum_{n=1}^N \exp(x_n - + * \max(x))\f$. + * + * @tparam T Type of input vector or matrix. + * @param[in] x Matrix of specified values. + * @return The log of the sum of the exponentiated vector values. + */ +template >>...> +inline auto log_sum_exp(T&& x) { + return apply_vector_unary::reduce(std::forward(x), [&](auto& v) { + using T_fvar_inner = typename value_type_t::Scalar; + using mat_type = Eigen::Matrix; + mat_type vals = v.val(); + mat_type exp_vals = vals.array().exp(); - return fvar(log_sum_exp(vals), - v.d().cwiseProduct(exp_vals).sum() / exp_vals.sum()); + return fvar( + log_sum_exp(vals), v.d().cwiseProduct(exp_vals).sum() / exp_vals.sum()); + }); } } // namespace math From 09c400415bede795037f753f63629e1742b89c8d Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Sun, 5 Jan 2020 10:38:32 +0000 Subject: [PATCH 15/22] [Jenkins] auto-formatting by clang-format version 5.0.0-3~16.04.1 (tags/RELEASE_500/final) --- stan/math/fwd/fun/log_sum_exp.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stan/math/fwd/fun/log_sum_exp.hpp b/stan/math/fwd/fun/log_sum_exp.hpp index 742045e594a..1881055c427 100644 --- a/stan/math/fwd/fun/log_sum_exp.hpp +++ b/stan/math/fwd/fun/log_sum_exp.hpp @@ -41,7 +41,7 @@ inline fvar log_sum_exp(const fvar& x1, double x2) { } template -fvar log_sum_exp(const std::vector >& v) { +fvar log_sum_exp(const std::vector>& v) { using std::exp; std::vector vals(v.size()); for (size_t i = 0; i < v.size(); ++i) { From a0eb3df5147379e8266dbb80a1a1713b94ebb8e3 Mon Sep 17 00:00:00 2001 From: rok-cesnovar Date: Sun, 5 Jan 2020 11:39:15 +0100 Subject: [PATCH 16/22] fix header guard --- stan/math/fwd/fun/log_softmax.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/stan/math/fwd/fun/log_softmax.hpp b/stan/math/fwd/fun/log_softmax.hpp index 7ba2ec759fb..2aa2fef9685 100644 --- a/stan/math/fwd/fun/log_softmax.hpp +++ b/stan/math/fwd/fun/log_softmax.hpp @@ -1,5 +1,5 @@ -#ifndef STAN_MATH_FWD_LOG_SOFTMAX_HPP -#define STAN_MATH_FWD_LOG_SOFTMAX_HPP +#ifndef STAN_MATH_FWD_FUN_LOG_SOFTMAX_HPP +#define STAN_MATH_FWD_FUN_LOG_SOFTMAX_HPP #include #include From 4e83afb8b183693f8f790c3b189879bb4f5e3302 Mon Sep 17 00:00:00 2001 From: rok-cesnovar Date: Sun, 5 Jan 2020 11:47:31 +0100 Subject: [PATCH 17/22] remove include --- stan/math/fwd/fun/log_sum_exp.hpp | 1 - 1 file changed, 1 deletion(-) diff --git a/stan/math/fwd/fun/log_sum_exp.hpp b/stan/math/fwd/fun/log_sum_exp.hpp index 1881055c427..63998a7fcd2 100644 --- a/stan/math/fwd/fun/log_sum_exp.hpp +++ b/stan/math/fwd/fun/log_sum_exp.hpp @@ -3,7 +3,6 @@ #include #include -#include #include #include #include From febffbeaa2d9fa91ead84b9ac9be1c9512a2d3c2 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Sat, 11 Jan 2020 12:49:15 +0800 Subject: [PATCH 18/22] Address review comments --- stan/math/fwd/fun/log_softmax.hpp | 4 +- stan/math/fwd/fun/log_sum_exp.hpp | 27 +------- stan/math/prim/mat/fun/log_softmax.hpp | 4 +- .../prim/vectorize/apply_vector_unary.hpp | 10 +-- stan/math/rev/core/matrix_vari.hpp | 6 +- stan/math/rev/fun/log_softmax.hpp | 4 +- stan/math/rev/fun/log_sum_exp.hpp | 61 ++----------------- test/unit/math/mix/fun/log_softmax_test.cpp | 28 ++++----- test/unit/math/mix/fun/log_sum_exp_test.cpp | 21 ++++--- 9 files changed, 48 insertions(+), 117 deletions(-) diff --git a/stan/math/fwd/fun/log_softmax.hpp b/stan/math/fwd/fun/log_softmax.hpp index 2aa2fef9685..1062038d2ad 100644 --- a/stan/math/fwd/fun/log_softmax.hpp +++ b/stan/math/fwd/fun/log_softmax.hpp @@ -20,8 +20,8 @@ namespace math { * @throw std::domain_error If the input vector is size 0. */ template >>...> -inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [&](auto& alpha) { +inline auto log_softmax(const T& x) { + return apply_vector_unary::apply(x, [&](auto& alpha) { using T_fvar = value_type_t; using T_fvar_inner = typename T_fvar::Scalar; diff --git a/stan/math/fwd/fun/log_sum_exp.hpp b/stan/math/fwd/fun/log_sum_exp.hpp index 63998a7fcd2..dcd0f401621 100644 --- a/stan/math/fwd/fun/log_sum_exp.hpp +++ b/stan/math/fwd/fun/log_sum_exp.hpp @@ -32,28 +32,7 @@ inline fvar log_sum_exp(double x1, const fvar& x2) { template inline fvar log_sum_exp(const fvar& x1, double x2) { - using std::exp; - if (x2 == NEGATIVE_INFTY) { - return fvar(x1.val_, x1.d_); - } - return fvar(log_sum_exp(x1.val_, x2), x1.d_ / (1 + exp(x2 - x1.val_))); -} - -template -fvar log_sum_exp(const std::vector>& v) { - using std::exp; - std::vector vals(v.size()); - for (size_t i = 0; i < v.size(); ++i) { - vals[i] = v[i].val_; - } - T deriv(0.0); - T denominator(0.0); - for (size_t i = 0; i < v.size(); ++i) { - T exp_vi = exp(vals[i]); - denominator += exp_vi; - deriv += v[i].d_ * exp_vi; - } - return fvar(log_sum_exp(vals), deriv / denominator); + return log_sum_exp(x2, x1); } /** @@ -72,8 +51,8 @@ fvar log_sum_exp(const std::vector>& v) { * @return The log of the sum of the exponentiated vector values. */ template >>...> -inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [&](auto& v) { +inline auto log_sum_exp(const T& x) { + return apply_vector_unary::reduce(x, [&](auto& v) { using T_fvar_inner = typename value_type_t::Scalar; using mat_type = Eigen::Matrix; mat_type vals = v.val(); diff --git a/stan/math/prim/mat/fun/log_softmax.hpp b/stan/math/prim/mat/fun/log_softmax.hpp index 6a6d8e5590e..4cd9dfb25fe 100644 --- a/stan/math/prim/mat/fun/log_softmax.hpp +++ b/stan/math/prim/mat/fun/log_softmax.hpp @@ -38,8 +38,8 @@ namespace math { * @return log unit simplex result of the softmax transform of the vector. */ template >>...> -inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [&](auto& v) { +inline auto log_softmax(const T& x) { + return apply_vector_unary::apply(x, [&](auto& v) { check_nonzero_size("log_softmax", "v", v); return (v.array() - log_sum_exp(v)).matrix(); }); diff --git a/stan/math/prim/vectorize/apply_vector_unary.hpp b/stan/math/prim/vectorize/apply_vector_unary.hpp index 9e5ad98386f..930969bb58d 100644 --- a/stan/math/prim/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/vectorize/apply_vector_unary.hpp @@ -45,7 +45,8 @@ struct apply_vector_unary> { /** * Member function for applying a functor to a vector and subsequently - * returning a scalar. + * returning a scalar. The reduction to a scalar needs to be implemented + * in the definition of the functor. * * @tparam T Type of argument to which functor is applied. * @tparam F Type of functor to apply. @@ -81,9 +82,10 @@ struct apply_vector_unary> { */ template static inline std::vector apply(const T& x, const F& f) { - Eigen::Matrix result - = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); - return std::vector(result.data(), result.data() + result.size()); + std::vector result(x.size()); + Eigen::Map>(result.data(), result.size()) + = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); + return result; } /** diff --git a/stan/math/rev/core/matrix_vari.hpp b/stan/math/rev/core/matrix_vari.hpp index 541ec5c7072..fe6d251c418 100644 --- a/stan/math/rev/core/matrix_vari.hpp +++ b/stan/math/rev/core/matrix_vari.hpp @@ -6,7 +6,6 @@ #include #include #include -#include namespace stan { namespace math { @@ -17,9 +16,8 @@ class op_matrix_vari : public vari { vari** vis_; public: - template >...> - op_matrix_vari(double f, const Eigen::MatrixBase& vs) - : vari(f), size_(vs.size()) { + template ...> + op_matrix_vari(double f, T&& vs) : vari(f), size_(vs.size()) { vis_ = ChainableStack::instance_->memalloc_.alloc_array(size_); Eigen::Map>(vis_, vs.rows(), vs.cols()) = vs.vi(); diff --git a/stan/math/rev/fun/log_softmax.hpp b/stan/math/rev/fun/log_softmax.hpp index 9ef071c0be1..c847ba848fc 100644 --- a/stan/math/rev/fun/log_softmax.hpp +++ b/stan/math/rev/fun/log_softmax.hpp @@ -55,8 +55,8 @@ class log_softmax_elt_vari : public vari { * @throw std::domain_error If the input vector is size 0. */ template >>...> -inline auto log_softmax(T&& x) { - return apply_vector_unary::apply(std::forward(x), [&](auto& alpha) { +inline auto log_softmax(const T& x) { + return apply_vector_unary::apply(x, [&](auto& alpha) { const int a_size = alpha.size(); check_nonzero_size("log_softmax", "alpha", alpha); diff --git a/stan/math/rev/fun/log_sum_exp.hpp b/stan/math/rev/fun/log_sum_exp.hpp index 30ea043cc21..d83eae9bbfb 100644 --- a/stan/math/rev/fun/log_sum_exp.hpp +++ b/stan/math/rev/fun/log_sum_exp.hpp @@ -5,10 +5,10 @@ #include #include #include -#include #include #include #include +#include #include #include @@ -38,18 +38,6 @@ class log_sum_exp_vd_vari : public op_vd_vari { } } }; -class log_sum_exp_dv_vari : public op_dv_vari { - public: - log_sum_exp_dv_vari(double a, vari* bvi) - : op_dv_vari(log_sum_exp(a, bvi->val_), a, bvi) {} - void chain() { - if (val_ == NEGATIVE_INFTY) { - bvi_->adj_ += adj_; - } else { - bvi_->adj_ += adj_ * calculate_chain(bvi_->val_, val_); - } - } -}; } // namespace internal @@ -69,45 +57,7 @@ inline var log_sum_exp(const var& a, double b) { * Returns the log sum of exponentials. */ inline var log_sum_exp(double a, const var& b) { - return var(new internal::log_sum_exp_dv_vari(a, b.vi_)); -} - -namespace internal { -inline double log_sum_exp_as_double(const std::vector& x) { - using std::exp; - using std::log; - double max = NEGATIVE_INFTY; - for (size_t i = 0; i < x.size(); ++i) { - if (x[i] > max) { - max = x[i].val(); - } - } - double sum = 0.0; - for (size_t i = 0; i < x.size(); ++i) { - if (x[i] != NEGATIVE_INFTY) { - sum += exp(x[i].val() - max); - } - } - return max + log(sum); -} - -class log_sum_exp_vector_vari : public op_vector_vari { - public: - explicit log_sum_exp_vector_vari(const std::vector& x) - : op_vector_vari(log_sum_exp_as_double(x), x) {} - void chain() { - for (size_t i = 0; i < size_; ++i) { - vis_[i]->adj_ += adj_ * calculate_chain(vis_[i]->val_, val_); - } - } -}; -} // namespace internal - -/** - * Returns the log sum of exponentials. - */ -inline var log_sum_exp(const std::vector& x) { - return var(new internal::log_sum_exp_vector_vari(x)); + return var(new internal::log_sum_exp_vd_vari(b.vi_, a)); } namespace internal { @@ -116,8 +66,7 @@ class log_sum_exp_matrix_vari : public op_matrix_vari { public: template explicit log_sum_exp_matrix_vari(T&& x) - : op_matrix_vari(log_sum_exp(std::forward(x).val()), - std::forward(x)) {} + : op_matrix_vari(log_sum_exp(x.val()), std::forward(x)) {} void chain() { Eigen::Map vis_map(vis_, size_); vis_map.adj().array() += adj_ * (vis_map.val().array() - val_).exp(); @@ -132,8 +81,8 @@ class log_sum_exp_matrix_vari : public op_matrix_vari { * @param x matrix */ template >>...> -inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [&](auto& v) { +inline auto log_sum_exp(const T& x) { + return apply_vector_unary::reduce(x, [&](auto& v) { return var(new internal::log_sum_exp_matrix_vari(v)); }); } diff --git a/test/unit/math/mix/fun/log_softmax_test.cpp b/test/unit/math/mix/fun/log_softmax_test.cpp index bee6b4ed61e..8155f9677fb 100644 --- a/test/unit/math/mix/fun/log_softmax_test.cpp +++ b/test/unit/math/mix/fun/log_softmax_test.cpp @@ -2,7 +2,7 @@ #include TEST(MathMixMatFun, logSoftmax) { - auto f = [](const auto& x) { return stan::math::log_softmax(x).eval(); }; + auto f = [](const auto& x) { return stan::math::log_softmax(x); }; // Column Vectors Eigen::VectorXd x0(0); // error case stan::test::expect_ad(f, x0); @@ -51,43 +51,41 @@ TEST(MathMixMatFun, logSoftmax) { rx3c << 2, 1, 1; stan::test::expect_ad(f, rx3c); - auto g = [](const auto& x) { return stan::math::log_softmax(x); }; - // std vectors std::vector stx0(0); // error case - stan::test::expect_ad(g, stx0); + stan::test::expect_ad(f, stx0); std::vector stx1{0}; - stan::test::expect_ad(g, stx1); + stan::test::expect_ad(f, stx1); std::vector stx2{-1, 1}; - stan::test::expect_ad(g, stx2); + stan::test::expect_ad(f, stx2); std::vector stx3{-1, 1, 10}; - stan::test::expect_ad(g, stx3); + stan::test::expect_ad(f, stx3); std::vector stx3b{0, 1, 2}; - stan::test::expect_ad(g, stx3b); + stan::test::expect_ad(f, stx3b); std::vector stx3c{2, 1, 1}; - stan::test::expect_ad(g, stx3c); + stan::test::expect_ad(f, stx3c); // Nested containers std::vector stvx0{x0, x0}; // error case - stan::test::expect_ad(g, stvx0); + stan::test::expect_ad(f, stvx0); std::vector stvx1{x1, x1}; - stan::test::expect_ad(g, stvx1); + stan::test::expect_ad(f, stvx1); std::vector strx0{rx0, rx0}; // error case - stan::test::expect_ad(g, strx0); + stan::test::expect_ad(f, strx0); std::vector strx1{rx1, rx1}; - stan::test::expect_ad(g, strx1); + stan::test::expect_ad(f, strx1); std::vector> ststx0{stx0, stx0}; // error case - stan::test::expect_ad(g, ststx0); + stan::test::expect_ad(f, ststx0); std::vector> ststx1{stx1, stx1}; - stan::test::expect_ad(g, ststx1); + stan::test::expect_ad(f, ststx1); } diff --git a/test/unit/math/mix/fun/log_sum_exp_test.cpp b/test/unit/math/mix/fun/log_sum_exp_test.cpp index aa3d9eaa46d..7244b70895e 100644 --- a/test/unit/math/mix/fun/log_sum_exp_test.cpp +++ b/test/unit/math/mix/fun/log_sum_exp_test.cpp @@ -70,18 +70,23 @@ TEST(MathMixMatFun, logSumExp) { stan::test::expect_ad(tols, f, x); Eigen::RowVectorXd rx = x; stan::test::expect_ad(tols, f, rx); + std::vector stx + = std::vector(x.data(), x.data() + x.size()); + stan::test::expect_ad(tols, f, stx); } Eigen::MatrixXd x23(2, 2); x23 << 1, 2, 3, 4; stan::test::expect_ad(f, x23); - std::vector a1{0}; - stan::test::expect_ad(f, a1); - - std::vector a2{5, 2}; - stan::test::expect_ad(f, a2); - - std::vector a4{1, 2, 3, 4}; - stan::test::expect_ad(f, a4); + std::vector stvx{x2, x2b, x2c}; + stan::test::expect_ad(tols, f, stvx); + std::vector strx{x2.transpose(), x2b.transpose(), + x2c.transpose()}; + stan::test::expect_ad(tols, f, strx); + std::vector> ststx{ + std::vector(x2.data(), x2.data() + x2.size()), + std::vector(x2b.data(), x2b.data() + x2b.size()), + std::vector(x2c.data(), x2c.data() + x2c.size())}; + stan::test::expect_ad(tols, f, ststx); } From 477cf9f12b96d57c1740c9b10820fa0d5d051e24 Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Sat, 11 Jan 2020 04:55:03 +0000 Subject: [PATCH 19/22] [Jenkins] auto-formatting by clang-format version 5.0.0-3~16.04.1 (tags/RELEASE_500/final) --- stan/math/prim/vectorize/apply_vector_unary.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stan/math/prim/vectorize/apply_vector_unary.hpp b/stan/math/prim/vectorize/apply_vector_unary.hpp index 930969bb58d..d2b95aa6105 100644 --- a/stan/math/prim/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/vectorize/apply_vector_unary.hpp @@ -84,7 +84,7 @@ struct apply_vector_unary> { static inline std::vector apply(const T& x, const F& f) { std::vector result(x.size()); Eigen::Map>(result.data(), result.size()) - = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); + = apply_vector_unary::apply(as_column_vector_or_scalar(x), f); return result; } From 02afdb9c98e047998c6d77cce5eb0aa26d00065d Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Sat, 11 Jan 2020 12:57:57 +0800 Subject: [PATCH 20/22] Fix merge error --- stan/math/prim/mat/fun/log_sum_exp.hpp | 4 +- .../math/mix/mat/fun/log_sum_exp_test.cpp | 54 ------------------- 2 files changed, 2 insertions(+), 56 deletions(-) delete mode 100644 test/unit/math/mix/mat/fun/log_sum_exp_test.cpp diff --git a/stan/math/prim/mat/fun/log_sum_exp.hpp b/stan/math/prim/mat/fun/log_sum_exp.hpp index fbaacf4a8e6..117d9814441 100644 --- a/stan/math/prim/mat/fun/log_sum_exp.hpp +++ b/stan/math/prim/mat/fun/log_sum_exp.hpp @@ -26,8 +26,8 @@ namespace math { * @return The log of the sum of the exponentiated vector values. */ template >>...> -inline auto log_sum_exp(T&& x) { - return apply_vector_unary::reduce(std::forward(x), [&](auto& v) { +inline auto log_sum_exp(const T& x) { + return apply_vector_unary::reduce(x, [&](auto& v) { if (v.size() == 0) { return NEGATIVE_INFTY; } diff --git a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp b/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp deleted file mode 100644 index c22d20a3aed..00000000000 --- a/test/unit/math/mix/mat/fun/log_sum_exp_test.cpp +++ /dev/null @@ -1,54 +0,0 @@ -#include -#include -#include - -TEST(MathMixMatFun, logSumExp) { - auto f = [](const auto& x) { return stan::math::log_sum_exp(x); }; - - Eigen::VectorXd x0(0); - stan::test::expect_ad(f, x0); - - Eigen::VectorXd x1(1); - x1 << 0; - - Eigen::VectorXd x2(2); - x2 << 5, 2; - - Eigen::VectorXd x2b(2); - x2b << 4.9, -std::numeric_limits::infinity(); - - Eigen::VectorXd x2c(2); - x2c << std::numeric_limits::infinity(), - std::numeric_limits::infinity(); - - Eigen::VectorXd x4(4); - x4 << 1, 2, 3, 4; - - stan::test::ad_tolerances tols; - tols.hessian_hessian_ = 1e-2; - tols.hessian_fvar_hessian_ = 1e-2; - - for (const auto& x : std::vector{x1, x2, x2b, x2c, x4}) { - stan::test::expect_ad(tols, f, x); - Eigen::RowVectorXd rx = x; - stan::test::expect_ad(tols, f, rx); - std::vector stx - = std::vector(x.data(), x.data() + x.size()); - stan::test::expect_ad(tols, f, stx); - } - - Eigen::MatrixXd x23(2, 2); - x23 << 1, 2, 3, 4; - stan::test::expect_ad(f, x23); - - std::vector stvx{x2, x2b, x2c}; - stan::test::expect_ad(tols, f, stvx); - std::vector strx{x2.transpose(), x2b.transpose(), - x2c.transpose()}; - stan::test::expect_ad(tols, f, strx); - std::vector> ststx{ - std::vector(x2.data(), x2.data() + x2.size()), - std::vector(x2b.data(), x2b.data() + x2b.size()), - std::vector(x2c.data(), x2c.data() + x2c.size())}; - stan::test::expect_ad(tols, f, ststx); -} From 8d8539ab09015314f9535bdab30da8c06c717599 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Sat, 11 Jan 2020 13:04:22 +0800 Subject: [PATCH 21/22] cpplint --- test/unit/math/prim/fun/log_sum_exp_test.cpp | 5 ----- 1 file changed, 5 deletions(-) diff --git a/test/unit/math/prim/fun/log_sum_exp_test.cpp b/test/unit/math/prim/fun/log_sum_exp_test.cpp index 67297bd073b..5876c8756a5 100644 --- a/test/unit/math/prim/fun/log_sum_exp_test.cpp +++ b/test/unit/math/prim/fun/log_sum_exp_test.cpp @@ -71,11 +71,6 @@ TEST(MathFunctions, log_sum_exp_nan) { EXPECT_TRUE(std::isnan(stan::math::log_sum_exp(nan, nan))); } -#include -#include -#include -#include - template void test_log_sum_exp(const Eigen::Matrix& as) { using stan::math::log_sum_exp; From f3b32861f683cc5a009027860e65ee26df966f32 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Mon, 13 Jan 2020 23:32:33 +0800 Subject: [PATCH 22/22] Address comments --- stan/math/fwd/fun/log_softmax.hpp | 2 +- stan/math/fwd/fun/log_sum_exp.hpp | 2 +- stan/math/prim/fun/log_softmax.hpp | 2 +- stan/math/prim/fun/log_sum_exp.hpp | 2 +- stan/math/prim/vectorize/apply_vector_unary.hpp | 2 +- stan/math/rev/core/matrix_vari.hpp | 2 +- stan/math/rev/fun/log_softmax.hpp | 2 +- stan/math/rev/fun/log_sum_exp.hpp | 6 +++--- 8 files changed, 10 insertions(+), 10 deletions(-) diff --git a/stan/math/fwd/fun/log_softmax.hpp b/stan/math/fwd/fun/log_softmax.hpp index 5782864faef..64d578cb4d7 100644 --- a/stan/math/fwd/fun/log_softmax.hpp +++ b/stan/math/fwd/fun/log_softmax.hpp @@ -22,7 +22,7 @@ namespace math { */ template >>...> inline auto log_softmax(const T& x) { - return apply_vector_unary::apply(x, [&](auto& alpha) { + return apply_vector_unary::apply(x, [&](const auto& alpha) { using T_fvar = value_type_t; using T_fvar_inner = typename T_fvar::Scalar; diff --git a/stan/math/fwd/fun/log_sum_exp.hpp b/stan/math/fwd/fun/log_sum_exp.hpp index 8d85462670b..1e8afb6f6e4 100644 --- a/stan/math/fwd/fun/log_sum_exp.hpp +++ b/stan/math/fwd/fun/log_sum_exp.hpp @@ -52,7 +52,7 @@ inline fvar log_sum_exp(const fvar& x1, double x2) { */ template >>...> inline auto log_sum_exp(const T& x) { - return apply_vector_unary::reduce(x, [&](auto& v) { + return apply_vector_unary::reduce(x, [&](const auto& v) { using T_fvar_inner = typename value_type_t::Scalar; using mat_type = Eigen::Matrix; mat_type vals = v.val(); diff --git a/stan/math/prim/fun/log_softmax.hpp b/stan/math/prim/fun/log_softmax.hpp index 21a7e03b8fd..9dde2af1f6d 100644 --- a/stan/math/prim/fun/log_softmax.hpp +++ b/stan/math/prim/fun/log_softmax.hpp @@ -39,7 +39,7 @@ namespace math { */ template >>...> inline auto log_softmax(const T& x) { - return apply_vector_unary::apply(x, [&](auto& v) { + return apply_vector_unary::apply(x, [&](const auto& v) { check_nonzero_size("log_softmax", "v", v); return (v.array() - log_sum_exp(v)).matrix(); }); diff --git a/stan/math/prim/fun/log_sum_exp.hpp b/stan/math/prim/fun/log_sum_exp.hpp index f57ae14244a..ded15b564cd 100644 --- a/stan/math/prim/fun/log_sum_exp.hpp +++ b/stan/math/prim/fun/log_sum_exp.hpp @@ -78,7 +78,7 @@ inline return_type_t log_sum_exp(const T2& a, const T1& b) { */ template >>...> inline auto log_sum_exp(const T& x) { - return apply_vector_unary::reduce(x, [&](auto& v) { + return apply_vector_unary::reduce(x, [&](const auto& v) { if (v.size() == 0) { return NEGATIVE_INFTY; } diff --git a/stan/math/prim/vectorize/apply_vector_unary.hpp b/stan/math/prim/vectorize/apply_vector_unary.hpp index d2b95aa6105..d485d5201c8 100644 --- a/stan/math/prim/vectorize/apply_vector_unary.hpp +++ b/stan/math/prim/vectorize/apply_vector_unary.hpp @@ -1,7 +1,7 @@ #ifndef STAN_MATH_PRIM_VECTORIZE_APPLY_VECTOR_UNARY_HPP #define STAN_MATH_PRIM_VECTORIZE_APPLY_VECTOR_UNARY_HPP -#include +#include #include #include diff --git a/stan/math/rev/core/matrix_vari.hpp b/stan/math/rev/core/matrix_vari.hpp index 3f67f16c7dc..4477d565bee 100644 --- a/stan/math/rev/core/matrix_vari.hpp +++ b/stan/math/rev/core/matrix_vari.hpp @@ -17,7 +17,7 @@ class op_matrix_vari : public vari { public: template ...> - op_matrix_vari(double f, T&& vs) : vari(f), size_(vs.size()) { + op_matrix_vari(double f, const T& vs) : vari(f), size_(vs.size()) { vis_ = ChainableStack::instance_->memalloc_.alloc_array(size_); Eigen::Map>(vis_, vs.rows(), vs.cols()) = vs.vi(); diff --git a/stan/math/rev/fun/log_softmax.hpp b/stan/math/rev/fun/log_softmax.hpp index 3a3c1ddef4d..a903419f738 100644 --- a/stan/math/rev/fun/log_softmax.hpp +++ b/stan/math/rev/fun/log_softmax.hpp @@ -57,7 +57,7 @@ class log_softmax_elt_vari : public vari { */ template >>...> inline auto log_softmax(const T& x) { - return apply_vector_unary::apply(x, [&](auto& alpha) { + return apply_vector_unary::apply(x, [&](const auto& alpha) { const int a_size = alpha.size(); check_nonzero_size("log_softmax", "alpha", alpha); diff --git a/stan/math/rev/fun/log_sum_exp.hpp b/stan/math/rev/fun/log_sum_exp.hpp index c4e04825597..471462ff149 100644 --- a/stan/math/rev/fun/log_sum_exp.hpp +++ b/stan/math/rev/fun/log_sum_exp.hpp @@ -65,8 +65,8 @@ namespace internal { class log_sum_exp_matrix_vari : public op_matrix_vari { public: template - explicit log_sum_exp_matrix_vari(T&& x) - : op_matrix_vari(log_sum_exp(x.val()), std::forward(x)) {} + explicit log_sum_exp_matrix_vari(const T& x) + : op_matrix_vari(log_sum_exp(x.val()), x) {} void chain() { Eigen::Map vis_map(vis_, size_); vis_map.adj().array() += adj_ * (vis_map.val().array() - val_).exp(); @@ -82,7 +82,7 @@ class log_sum_exp_matrix_vari : public op_matrix_vari { */ template >>...> inline auto log_sum_exp(const T& x) { - return apply_vector_unary::reduce(x, [&](auto& v) { + return apply_vector_unary::reduce(x, [&](const auto& v) { return var(new internal::log_sum_exp_matrix_vari(v)); }); }