Get all the unit tests passing with the new API

stan/math/fwd/mat/meta/operands_and_partials.hpp

@@ -9,35 +9,34 @@
 namespace stan {
   namespace math {
     namespace detail {
-      template <typename ViewElt, typename Op, typename Dx>
-      class ops_partials_edge_vec_fwd
-        : public ops_partials_edge_mat_prim<ViewElt, Op, -1, 1> {
+      // Vectorized Univariate
+      template <typename ViewElt, typename Dx>
+      class ops_partials_edge<ViewElt, std::vector<fvar<Dx> > >
+        : public ops_partials_edge_mat_prim<ViewElt, std::vector<fvar<Dx> >, -1, 1> {
       public:
-        ops_partials_edge_vec_fwd(const Op& ops)
-          : ops_partials_edge_mat_prim<ViewElt, Op, -1, 1>(ops) {}
+        explicit ops_partials_edge(const std::vector<fvar<Dx> >& ops)
+          : ops_partials_edge_mat_prim<ViewElt, std::vector<fvar<Dx> >, -1, 1>(ops) {}
         Dx dx() {
           Dx derivative(0);
           for (int i = 0; i < this->size(); ++i) {
-            derivative += this->partials[i] * this->operands[i].d_;
+            derivative += this->partials(i) * this->operands[i].d_;
           }
           return derivative;
         }
       };
-
-      // Vectorized Univariate
-      template <typename ViewElt, typename Dx>
-      class ops_partials_edge<ViewElt, std::vector<fvar<Dx> > >
-        : public ops_partials_edge_vec_fwd<ViewElt, std::vector<fvar<Dx> >, Dx> {
-      public:
-        explicit ops_partials_edge(const std::vector<fvar<Dx> >& ops)
-          : ops_partials_edge_vec_fwd<ViewElt, std::vector<fvar<Dx> >, Dx>(ops) {}
-      };
       template <typename ViewElt, typename Dx, int R, int C>
       class ops_partials_edge<ViewElt, Eigen::Matrix<fvar<Dx>, R, C> >
-        : public ops_partials_edge_vec_fwd<ViewElt, Eigen::Matrix<fvar<Dx>, R, C>, Dx> {
+        : public ops_partials_edge_mat_prim<ViewElt, Eigen::Matrix<fvar<Dx>, R, C>, R, C> {
       public:
         explicit ops_partials_edge(const Eigen::Matrix<fvar<Dx>, R, C>& ops)
-          : ops_partials_edge_vec_fwd<ViewElt, Eigen::Matrix<fvar<Dx>, R, C>, Dx>(ops) {}
+          : ops_partials_edge_mat_prim<ViewElt, Eigen::Matrix<fvar<Dx>, R, C>, R, C>(ops) {}
+        Dx dx() {
+          Dx derivative(0);
+          for (int i = 0; i < this->size(); ++i) {
+            derivative += this->partials(i) * this->operands(i).d_;
+          }
+          return derivative;
+        }
       };
 
       // Multivariate; vectors of eigen types

stan/math/fwd/scal/meta/OperandsAndPartials.hpp

-#ifndef STAN_MATH_FWD_SCAL_META_OPERANDSANDPARTIALS_HPP
-#define STAN_MATH_FWD_SCAL_META_OPERANDSANDPARTIALS_HPP
-
-#include <stan/math/prim/scal/meta/is_constant_struct.hpp>
-#include <stan/math/prim/scal/meta/length.hpp>
-#include <stan/math/prim/scal/meta/OperandsAndPartials.hpp>
-#include <stan/math/fwd/core.hpp>
-
-namespace stan {
-  namespace math {
-
-    namespace {
-      template <typename T_derivative,
-                typename T,
-                typename T_partials,
-                bool is_vec = is_vector<T>::value,
-                bool is_const = is_constant_struct<T>::value>
-      struct increment_derivative {
-        inline T_derivative operator()(const T& x,
-                                       const T_partials& d_dx) {
-          return 0;
-        }
-      };
-
-      template <typename T_derivative,
-                typename T,
-                typename T_partials>
-      struct increment_derivative<T_derivative, T, T_partials, false, false> {
-        inline T_derivative operator()(const T& x,
-                                       const T_partials& d_dx) {
-          return d_dx[0] * x.d_;
-        }
-      };
-
-      template <typename T_derivative,
-                typename T,
-                typename T_partials>
-      struct increment_derivative<T_derivative, T, T_partials, true, false> {
-        inline T_derivative operator()(const T& x,
-                                       const T_partials& d_dx) {
-          T_derivative derivative(0);
-          for (size_t n = 0; n < length(x); n++)
-            derivative += d_dx[n] * x[n].d_;
-          return derivative;
-        }
-      };
-
-      template <typename T,
-                typename T1, typename D1, typename T2, typename D2,
-                typename T3, typename D3, typename T4, typename D4,
-                typename T5, typename D5, typename T6, typename D6>
-      fvar<T> partials_to_fvar(T& logp,
-                               const T1& x1, D1& d_x1,
-                               const T2& x2, D2& d_x2,
-                               const T3& x3, D3& d_x3,
-                               const T4& x4, D4& d_x4,
-                               const T5& x5, D5& d_x5,
-                               const T6& x6, D6& d_x6) {
-        T deriv = 0;
-        if (!is_constant_struct<T1>::value)
-          deriv += increment_derivative<T, T1, D1>()(x1, d_x1);
-        if (!is_constant_struct<T2>::value)
-          deriv += increment_derivative<T, T2, D2>()(x2, d_x2);
-        if (!is_constant_struct<T3>::value)
-          deriv += increment_derivative<T, T3, D3>()(x3, d_x3);
-        if (!is_constant_struct<T4>::value)
-          deriv += increment_derivative<T, T4, D4>()(x4, d_x4);
-        if (!is_constant_struct<T5>::value)
-          deriv += increment_derivative<T, T5, D5>()(x5, d_x5);
-        if (!is_constant_struct<T6>::value)
-          deriv += increment_derivative<T, T6, D6>()(x6, d_x6);
-        return fvar<T>(logp, deriv);
-      }
-    }
-
-    /**
-     * This class builds partial derivatives with respect to a set of
-     * operands. There are two reason for the generality of this
-     * class. The first is to handle vector and scalar arguments
-     * without needing to write additional code. The second is to use
-     * this class for writing probability distributions that handle
-     * primitives, reverse mode, and forward mode variables
-     * seamlessly.
-     *
-     * This is the partial template specialization for when the return
-     * type is fvar<T>.
-     *
-     * @tparam T1 First set of operands.
-     * @tparam T2 Second set of operands.
-     * @tparam T3 Third set of operands.
-     * @tparam T4 Fourth set of operands.
-     * @tparam T5 Fifth set of operands.
-     * @tparam T6 Sixth set of operands.
-     * @tparam T_return_type Return type of the expression. This defaults
-     *   to a template metaprogram that calculates the scalar promotion of
-     *   T1 -- T6.
-     */
-    template<typename T1, typename T2, typename T3,
-             typename T4, typename T5, typename T6,
-             typename T_partials_return>
-    struct OperandsAndPartials<T1, T2, T3, T4, T5, T6,
-                               fvar<T_partials_return> > {
-      typedef fvar<T_partials_return> T_return_type;
-
-      const T1& x1_;
-      const T2& x2_;
-      const T3& x3_;
-      const T4& x4_;
-      const T5& x5_;
-      const T6& x6_;
-
-      size_t n_partials;
-      T_partials_return* all_partials;
-
-      VectorView<T_partials_return,
-                 is_vector<T1>::value,
-                 is_constant_struct<T1>::value> d_x1;
-      VectorView<T_partials_return,
-                 is_vector<T2>::value,
-                 is_constant_struct<T2>::value> d_x2;
-      VectorView<T_partials_return,
-                 is_vector<T3>::value,
-                 is_constant_struct<T3>::value> d_x3;
-      VectorView<T_partials_return,
-                 is_vector<T4>::value,
-                 is_constant_struct<T4>::value> d_x4;
-      VectorView<T_partials_return,
-                 is_vector<T5>::value,
-                 is_constant_struct<T5>::value> d_x5;
-      VectorView<T_partials_return,
-                 is_vector<T6>::value,
-                 is_constant_struct<T6>::value> d_x6;
-
-      OperandsAndPartials(const T1& x1 = 0, const T2& x2 = 0, const T3& x3 = 0,
-                          const T4& x4 = 0, const T5& x5 = 0, const T6& x6 = 0)
-        : x1_(x1), x2_(x2), x3_(x3), x4_(x4), x5_(x5), x6_(x6),
-          n_partials(!is_constant_struct<T1>::value * length(x1) +
-                     !is_constant_struct<T2>::value * length(x2) +
-                     !is_constant_struct<T3>::value * length(x3) +
-                     !is_constant_struct<T4>::value * length(x4) +
-                     !is_constant_struct<T5>::value * length(x5) +
-                     !is_constant_struct<T6>::value * length(x6)),
-          all_partials(new T_partials_return[n_partials]),
-          d_x1(all_partials),
-          d_x2(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)),
-          d_x3(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)),
-          d_x4(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)
-               + (!is_constant_struct<T3>::value) * length(x3)),
-          d_x5(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)
-               + (!is_constant_struct<T3>::value) * length(x3)
-               + (!is_constant_struct<T4>::value) * length(x4)),
-          d_x6(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)
-               + (!is_constant_struct<T3>::value) * length(x3)
-               + (!is_constant_struct<T4>::value) * length(x4)
-               + (!is_constant_struct<T5>::value) * length(x5)) {
-        std::fill(all_partials, all_partials + n_partials, 0);
-      }
-
-      T_return_type
-      value(T_partials_return value) {
-        return partials_to_fvar(value,
-                                x1_, d_x1, x2_, d_x2,
-                                x3_, d_x3, x4_, d_x4,
-                                x5_, d_x4, x6_, d_x5);
-      }
-
-      ~OperandsAndPartials() {
-        delete[] all_partials;
-      }
-    };
-
-  }
-}
-#endif

stan/math/fwd/scal/meta/operands_and_partials.hpp

@@ -15,80 +15,40 @@ namespace stan {
         ops_partials_edge(const fvar<Dx>& op)
           : ops_partials_edge_singular<ViewElt, fvar<Dx> >(op) {}
         Dx dx() {
-          return this->partial * this->operand.d_;
+          return this->partials[0] * this->operand.d_;
         }
       };
+    } // end namespace detail
+    template <typename Op1, typename Op2, typename Op3, typename Op4,
+              typename Dx>
+    class operands_and_partials<Op1, Op2, Op3, Op4, fvar<Dx> > {
+    public:
+      // these are going to be stack local and get collected
+      detail::ops_partials_edge<Dx, Op1> edge1_;
+      detail::ops_partials_edge<Dx, Op2> edge2_;
+      detail::ops_partials_edge<Dx, Op3> edge3_;
+      detail::ops_partials_edge<Dx, Op4> edge4_;
+      typedef fvar<Dx> T_return_type;
+      operands_and_partials(const Op1& o1)
+        : edge1_(o1) { }
+      operands_and_partials(const Op1& o1, const Op2& o2)
+        : edge1_(o1), edge2_(o2) { }
+      operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3)
+        : edge1_(o1), edge2_(o2), edge3_(o3) { }
+      operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3, const Op4& o4)
+        : edge1_(o1), edge2_(o2), edge3_(o3), edge4_(o4) { }
 
-      template <typename Op1, typename Op2, typename Op3, typename Op4,
-                typename Dx>
-      class operands_and_partials<Op1, Op2, Op3, Op4, fvar<Dx> > {
-      public:
-        typedef fvar<Dx> T_return_type;
-        operands_and_partials(const Op1& o1)
-          : edge1_(o1) { }
-        operands_and_partials(const Op1& o1, const Op2& o2)
-          : edge1_(o1), edge2_(o2) { }
-        operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3)
-          : edge1_(o1), edge2_(o2), edge3_(o3) { }
-        operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3, const Op4& o4)
-          : edge1_(o1), edge2_(o2), edge3_(o3), edge4_(o4) { }
-
-        void increment_dx1(const int n, const Dx& adj) {
-          edge1_.increment_dx(n, adj);
-        }
-        void increment_dx2(const int n, const Dx& adj) {
-          edge2_.increment_dx(n, adj);
-        }
-        void increment_dx3(const int n, const Dx& adj) {
-          edge3_.increment_dx(n, adj);
-        }
-        void increment_dx4(const int n, const Dx& adj) {
-          edge4_.increment_dx(n, adj);
-        }
-
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op1>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx1_vector(const int n, const T& adj) {
-          edge1_.increment_dx_vector(n, adj);
-        }
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op2>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx2_vector(const int n, const T& adj) {
-          edge2_.increment_dx_vector(n, adj);
-        }
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op3>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx3_vector(const int n, const T& adj) {
-          edge3_.increment_dx_vector(n, adj);
-        }
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op4>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx4_vector(const int n, const T& adj) {
-          edge4_.increment_dx_vector(n, adj);
-        }
-
-        // this is what matters in terms of going on autodiff stack
-        T_return_type build(Dx value) {
-          // Just call chain here and end up with an fvar(val, adj)
-          Dx deriv = 0;
-          deriv += edge1_.dx();
-          deriv += edge2_.dx();
-          deriv += edge3_.dx();
-          deriv += edge4_.dx();
-          return T_return_type(value, deriv);
-        };
-      private:
-        // these are going to be stack local and get collected
-        ops_partials_edge<Dx, Op1> edge1_;
-        ops_partials_edge<Dx, Op2> edge2_;
-        ops_partials_edge<Dx, Op3> edge3_;
-        ops_partials_edge<Dx, Op4> edge4_;
+      // this is what matters in terms of going on autodiff stack
+      T_return_type build(Dx value) {
+        // Just call chain here and end up with an fvar(val, adj)
+        Dx deriv = 0;
+        deriv += edge1_.dx();
+        deriv += edge2_.dx();
+        deriv += edge3_.dx();
+        deriv += edge4_.dx();
+        return T_return_type(value, deriv);
       };
-    } // end namespace detail
+    };
   }
 }
 

stan/math/prim/arr/meta/VectorView.hpp

-#ifndef STAN_MATH_ARR_SCAL_META_VECTORVIEW_HPP
-#define STAN_MATH_ARR_SCAL_META_VECTORVIEW_HPP
-
-#include <stan/math/prim/scal/meta/VectorView.hpp>
-#include <vector>
-
-namespace stan {
-
-  template <typename T>
-  class VectorView<std::vector<T>, true, false> {
-  public:
-    typedef typename scalar_type<T>::type scalar_t;
-
-    template <typename X>
-    explicit VectorView(X& x) : x_(&x[0]) { }
-
-    scalar_t& operator[](int i) {
-      return x_[i];
-    }
-
-  private:
-    scalar_t* x_;
-  };
-
-  template <typename T>
-  class VectorView<const std::vector<T>, true, false> {
-  public:
-    typedef typename boost::add_const<typename scalar_type<T>::type>::type
-    scalar_t;
-
-    template <typename X>
-    explicit VectorView(X& x) : x_(&x[0]) { }
-
-    scalar_t& operator[](int i) const {
-      return x_[i];
-    }
-  private:
-    scalar_t* x_;
-  };
-
-}
-#endif

stan/math/prim/mat/meta/VectorView.hpp

-#ifndef STAN_MATH_MAT_SCAL_META_VECTORVIEW_HPP
-#define STAN_MATH_MAT_SCAL_META_VECTORVIEW_HPP
-
-#include <stan/math/prim/mat/fun/Eigen.hpp>
-#include <stan/math/prim/mat/meta/scalar_type.hpp>
-#include <stan/math/prim/scal/meta/VectorView.hpp>
-#include <boost/type_traits.hpp>
-
-namespace stan {
-
-  template <typename T, int R, int C>
-  class VectorView<Eigen::Matrix<T, R, C>, true, false> {
-  public:
-    typedef typename scalar_type<T>::type scalar_t;
-
-    template <typename X>
-    explicit VectorView(X& x) : x_(x.data()) { }
-
-    scalar_t& operator[](int i) {
-      return x_[i];
-    }
-  private:
-    scalar_t* x_;
-  };
-
-  template <typename T, int R, int C>
-  class VectorView<const Eigen::Matrix<T, R, C>, true, false> {
-  public:
-    typedef typename boost::add_const<typename scalar_type<T>::type>::type
-    scalar_t;
-
-    template <typename X>
-    explicit VectorView(X& x) : x_(x.data()) { }
-
-    scalar_t& operator[](int i) const {
-      return x_[i];
-    }
-  private:
-    scalar_t* x_;
-  };
-
-}
-#endif

stan/math/prim/mat/meta/operands_and_partials.hpp

@@ -3,30 +3,48 @@
 
 #include <stan/math/prim/mat/fun/typedefs.hpp>
 #include <stan/math/prim/scal/meta/operands_and_partials.hpp>
+#include <stan/math/prim/scal/meta/broadcast_array.hpp>
 #include <vector>
 
 namespace stan {
   namespace math {
     namespace detail {
+      template <typename T, typename Orig, int R, int C>
+      struct zero_vec_or_mat {
+        typedef Eigen::Matrix<T, R, C> ret;
+        static ret zero(const Orig& in) {
+          return ret::Zero(in.rows(), in.cols());
+        }
+      };
+      template <typename T, typename Orig>
+      struct zero_vec_or_mat<T, Orig, -1, 1> {
+        typedef Eigen::Matrix<T, -1, 1> ret;
+        static ret zero(const Orig& in) {
+          return ret::Zero(in.size(), 1);
+        }
+      };
+      template <typename T, typename Orig>
+      struct zero_vec_or_mat<T, Orig, 1, -1> {
+        typedef Eigen::Matrix<T, 1, -1> ret;
+        static ret zero(const Orig& in) {
+          return ret::Zero(1, in.size());
+        }
+      };
+
       // Handles all the vectorized cases for "views of scalars":
       // ViewElt = double, Arg = std::vector<var> d_ holds vector<double>
       // ViewElt = double, Arg = std::vector<double> d_ holds dummy
       // Op will always be some container of vars
       template <typename ViewElt, typename Op, int R, int C>
       class ops_partials_edge_mat_prim {
-      protected:
+      public:
         typedef Eigen::Matrix<ViewElt, R, C> partials_t;
         const Op& operands;
         partials_t partials;
-      public:
+        broadcast_array<partials_t> partials_vec;
         ops_partials_edge_mat_prim(const Op& ops)
-          : operands(ops), partials(partials_t::Zero(ops.size())) {}
-        void increment_dx(int n, const ViewElt& adj) {
-          partials[n] += adj;
-        }
-        void increment_dx_vector(int /*n*/, const partials_t& adj) {
-          partials += adj;
-        }
+          : operands(ops), partials(zero_vec_or_mat<ViewElt, Op, R, C>::zero(ops)),
+            partials_vec(partials) {}
         void dump_partials(double* partials) {
           for (int i = 0; i < this->partials.size(); ++i) {
             partials[i] = this->partials(i);
@@ -41,21 +59,17 @@ namespace stan {
       class ops_partials_edge_multivariate_prim {
       public:
         typedef Eigen::Matrix<ViewElt, Eigen::Dynamic, Eigen::Dynamic> partial_t;
+        std::vector<partial_t> partials_vec;
+        const std::vector<Op>& operands;
         ops_partials_edge_multivariate_prim(const std::vector<Op>& ops)
-          : partials(ops.size()), operands(ops) {
+          : partials_vec(ops.size()), operands(ops) {
           for (size_t i = 0; i < ops.size(); ++i) {
-            partials[i] = partial_t::Zero(ops[i].rows(), ops[i].cols());
+            partials_vec[i] = partial_t::Zero(ops[i].rows(), ops[i].cols());
           }
         }
-        void increment_dx_vector(int n, const partial_t& adj) {
-          partials[n] += adj;
-        }
         int size() {
           return this->operands.size() * this->operands[0].size();
         }
-      protected:
-        std::vector<partial_t> partials;
-        const std::vector<Op>& operands;
       };
     }
   }

stan/math/prim/scal/meta/OperandsAndPartials.hpp

-#ifndef STAN_MATH_PRIM_SCAL_META_OPERANDSANDPARTIALS_HPP
-#define STAN_MATH_PRIM_SCAL_META_OPERANDSANDPARTIALS_HPP
-
-#include <stan/math/prim/scal/meta/return_type.hpp>
-#include <stan/math/prim/scal/meta/VectorView.hpp>
-
-namespace stan {
-  namespace math {
-
-    /**
-     * This class builds partial derivatives with respect to a set of
-     * operands. There are two reason for the generality of this
-     * class. The first is to handle vector and scalar arguments
-     * without needing to write additional code. The second is to use
-     * this class for writing probability distributions that handle
-     * primitives, reverse mode, and forward mode variables
-     * seamlessly.
-     * 
-     * The default template class handles the case where the arguments
-     * are primitive. There are template specializations for reverse
-     * mode and forward mode.
-     *
-     * @tparam T1 First set of operands.
-     * @tparam T2 Second set of operands.
-     * @tparam T3 Third set of operands.
-     * @tparam T4 Fourth set of operands.
-     * @tparam T5 Fifth set of operands.
-     * @tparam T6 Sixth set of operands.
-     * @tparam T_return_type Return type of the expression. This defaults
-     *   to a template metaprogram that calculates the scalar promotion of
-     *   T1 -- T6.
-     */
-    template<typename T1 = double, typename T2 = double, typename T3 = double,
-             typename T4 = double, typename T5 = double, typename T6 = double,
-             typename T_return_type
-             = typename stan::return_type<T1, T2, T3, T4, T5, T6>::type>
-    struct OperandsAndPartials {
-      VectorView<T_return_type, false, true> d_x1;
-      VectorView<T_return_type, false, true> d_x2;
-      VectorView<T_return_type, false, true> d_x3;
-      VectorView<T_return_type, false, true> d_x4;
-      VectorView<T_return_type, false, true> d_x5;
-      VectorView<T_return_type, false, true> d_x6;
-
-      /**
-       * Constructor.
-       *
-       * @param x1 first set of operands
-       * @param x2 second set of operands
-       * @param x3 third set of operands
-       * @param x4 fourth set of operands
-       * @param x5 fifth set of operands
-       * @param x6 sixth set of operands
-       */
-      OperandsAndPartials(const T1& x1 = 0, const T2& x2 = 0,
-                          const T3& x3 = 0, const T4& x4 = 0,
-                          const T5& x5 = 0, const T6& x6 = 0) { }
-
-      /**
-       * Returns a T_return_type with the value specified with
-       * the partial derivatves.
-       *
-       * @param[in] value Value of the variable
-       * @returns a variable with the appropriate value
-       */
-      T_return_type
-      value(double value) {
-        return value;
-      }
-    };
-
-  }
-}
-#endif

stan/math/prim/scal/meta/VectorView.hpp

-#ifndef STAN_MATH_PRIM_SCAL_META_VECTORVIEW_HPP
-#define STAN_MATH_PRIM_SCAL_META_VECTORVIEW_HPP
-
-#include <stan/math/prim/scal/meta/scalar_type.hpp>
-#include <stan/math/prim/scal/meta/is_vector_like.hpp>
-#include <boost/type_traits.hpp>
-#include <stdexcept>
-
-namespace stan {
-
-  /**
-   * VectorView is a template expression that is constructed with a
-   * container or scalar, which it then allows to be used as an array
-   * using <code>operator[]</code>.
-   *
-   * For a scalar value, any index returns the reference or pointer
-   * used to construct the view.
-   *
-   * For a container, the index returns a reference to the position in
-   * the underlying container used to construct the view.  WARNING:
-   * There is no bounds checking for container indices and they will
-   * segfault if accessed beyond their boundaries.
-   *
-   * The first use is to read arguments to prob functions as vectors,
-   * even if scalars, so they can be read by common code (and scalars
-   * automatically broadcast up to behave like vectors) : VectorView
-   * of immutable const array of double* (no allocation).
-   *
-   * The second use is to build up derivatives into common storage :
-   * VectorView of mutable shared array (no allocation because
-   * allocated on auto-diff arena memory).
-   *
-   * Because it deals with references to its inputs, it is up to the
-   * client of VectorView to ensure that the container being wrapped
-   * is not modified while the VectorView is in use in such a way as
-   * to disrupt the indexing.  Similarly, because it deals with
-   * references, it cannot be constructed with a literal or expression.
-   *
-   * @tparam T  Type of scalar or container being wrapped.
-   * @tparam is_array True if underlying type T can be indexed with
-   * operator[].
-   * @tparam throw_if_accessed True if the behavior is to throw an
-   * exception whenever <code>operator[]</code> is called.
-   */
-  template <typename T,
-            bool is_array = stan::is_vector_like<T>::value,
-            bool throw_if_accessed = false>
-  class VectorView {
-  public:
-    typedef typename
-    boost::conditional<boost::is_const<T>::value,
-                       typename boost::add_const<
-                         typename scalar_type<T>::type>::type,
-                       typename scalar_type<T>::type>::type scalar_t;
-
-    template <typename X>
-    explicit VectorView(X x) {
-      throw std::logic_error("VectorView: the default template "
-                             "specialization not implemented");
-    }
-
-    scalar_t& operator[](int i) {
-      throw std::logic_error("VectorView: the default template "
-                             "specialization not implemented");
-    }
-
-    scalar_t& operator[](int i) const {
-      throw std::logic_error("VectorView: the default template "
-                             "specialization not implemented");
-    }
-  };
-
-  template <typename T, bool is_array>
-  class VectorView<T, is_array, true> {
-  public:
-    typedef typename
-    boost::conditional<boost::is_const<T>::value,
-                       typename boost::add_const<
-                         typename scalar_type<T>::type>::type,
-                       typename scalar_type<T>::type>::type scalar_t;
-    VectorView() { }
-
-    template <typename X>
-    explicit VectorView(X x) { }
-
-    scalar_t& operator[](int i) {
-      throw std::logic_error("VectorView: this cannot be accessed");
-    }
-
-    scalar_t& operator[](int i) const {
-      throw std::logic_error("VectorView: this cannot be accessed");
-    }
-  };
-
-  // this covers non-vectors: double
-  template <typename T>
-  class VectorView<T, false, false> {
-  public:
-    typedef typename
-    boost::conditional<boost::is_const<T>::value,
-                       typename boost::add_const<
-                         typename scalar_type<T>::type>::type,
-                       typename scalar_type<T>::type>::type scalar_t;
-
-    explicit VectorView(scalar_t& x) : x_(&x) { }
-
-    explicit VectorView(scalar_t* x) : x_(x) { }
-
-    scalar_t& operator[](int i) {
-      return *x_;
-    }
-
-    scalar_t& operator[](int i) const {
-      return *x_;
-    }
-  private:
-    scalar_t* x_;
-  };
-
-  // this covers raw memory: double*
-  template <typename T>
-  class VectorView<T, true, false> {
-  public:
-    typedef typename
-    boost::conditional<boost::is_const<T>::value,
-                       typename boost::add_const<
-                         typename scalar_type<T>::type>::type,
-                       typename scalar_type<T>::type>::type scalar_t;
-
-    explicit VectorView(scalar_t* x) : x_(x) { }
-
-    scalar_t& operator[](int i) {
-      return x_[i];
-    }
-
-    scalar_t& operator[](int i) const {
-      return x_[i];
-    }
-
-  private:
-    scalar_t* x_;
-  };
-}
-#endif

stan/math/prim/scal/meta/broadcast_array.hpp

+#ifndef STAN_MATH_PRIM_SCAL_META_BROADCAST_ARRAY_HPP
+#define STAN_MATH_PRIM_SCAL_META_BROADCAST_ARRAY_HPP
+
+namespace stan {
+  namespace math {
+    namespace detail {
+      template <typename T>
+      class broadcast_array {
+      private:
+        T& prim_;
+
+      public:
+        explicit broadcast_array(T& prim): prim_(prim) {}
+
+        T& operator[] (int /*i*/) {
+          return prim_;
+        }
+        int size() { return 1; }
+      };
+
+      template <>
+      class broadcast_array<void> {
+      public:
+        explicit broadcast_array() {}
+        double& operator[] (int /*i*/) {
+          throw std::logic_error("Don't do this");
+        }
+        int size() { return 0; }
+      };
+    } // end namespace detail
+  } // end namespace math
+} // end namespace stan
+
+#endif

stan/math/prim/scal/meta/operands_and_partials.hpp

@@ -2,6 +2,7 @@
 #define STAN_MATH_PRIM_SCAL_META_OPERANDS_AND_PARTIALS_HPP
 
 #include <stan/math/prim/scal/meta/return_type.hpp>
+#include <stan/math/prim/scal/meta/broadcast_array.hpp>
 
 namespace stan {
   namespace math {
@@ -52,40 +53,10 @@ namespace stan {
       public:
         ops_partials_edge() { }
         ops_partials_edge(const Op& /* a */){ }
-        void increment_dx(int /* n */, const ViewElt& /* adj */) {}
         void dump_partials(ViewElt* /* partials */) {} // used for vars
         void dump_operands(void* /* operands */) {} // also used for reverse mode
         double dx() const { return 0; } //used for fvars
-        int size() {return 0; }
-      };
-
-      template <typename Op1 = double, typename Op2 = double,
-                typename Op3 = double, typename Op4 = double,
-                typename T_return_type = typename return_type<Op1, Op2, Op3, Op4>::type>
-      class operands_and_partials {
-      public:
-        operands_and_partials(const Op1& op1) {}
-        operands_and_partials(const Op1& op1, const Op2& op2) {}
-        operands_and_partials(const Op1& op1, const Op2& op2, const Op3& op3) {}
-        operands_and_partials(const Op1& op1, const Op2& op2, const Op3& op3, const Op4& op4) {}
-
-        void increment_dx1(int, double) {}
-        void increment_dx2(int, double) {}
-        void increment_dx3(int, double) {}
-        void increment_dx4(int, double) {}
-
-        template <typename PartialVec>
-        void increment_dx1_vector(int, PartialVec) {}
-        template <typename PartialVec>
-        void increment_dx2_vector(int, PartialVec) {}
-        template <typename PartialVec>
-        void increment_dx3_vector(int, PartialVec) {}
-        template <typename PartialVec>
-        void increment_dx4_vector(int, PartialVec) {}
-
-        double build(const double value) {
-          return value;
-        }
+        int size() { return 0; }
       };
 
       // Base class shared between fvar and var specializations of uni- and
@@ -93,21 +64,45 @@ namespace stan {
       // Singular version - underlying Var is the same as Op, so there's just one.
       template <typename ViewElt, typename Op>
       class ops_partials_edge_singular {
-      protected:
+      public:
         const Op& operand;
         ViewElt partial;
-
-      public:
-        ops_partials_edge_singular(const Op& a) : operand(a), partial(0) { }
-        void increment_dx(int /*n*/, const ViewElt& adj) {
-          partial += adj;
-        }
+        broadcast_array<ViewElt> partials;
+        ops_partials_edge_singular(const Op& a) : operand(a), partial(0), partials(partial) { }
         // dump_operands implemented in specialization
-        int size() {
-          return 1;
-        }
+        int size() { return 1; }
+      };
+
+      template<>
+      struct ops_partials_edge<double, double> {
+        ops_partials_edge() {}
+        ops_partials_edge(const double& /* a */) {}
+        broadcast_array<void> partials;
+        void dump_operands(void*) {}
+        void dump_partials(void*) {}
+        double dx() { return 0; }
+        int size() { return 0; }
       };
     } // end namespace detail
+    template <typename Op1 = double, typename Op2 = double,
+              typename Op3 = double, typename Op4 = double,
+              typename T_return_type = typename return_type<Op1, Op2, Op3, Op4>::type>
+    class operands_and_partials {
+    public:
+      operands_and_partials(const Op1& op1) {}
+      operands_and_partials(const Op1& op1, const Op2& op2) {}
+      operands_and_partials(const Op1& op1, const Op2& op2, const Op3& op3) {}
+      operands_and_partials(const Op1& op1, const Op2& op2, const Op3& op3, const Op4& op4) {}
+
+      double build(const double value) {
+        return value;
+      }
+
+      detail::ops_partials_edge<double, Op1> edge1_;
+      detail::ops_partials_edge<double, Op2> edge2_;
+      detail::ops_partials_edge<double, Op3> edge3_;
+      detail::ops_partials_edge<double, Op4> edge4_;
+    };
   } // end namespace math
 } // end namespace stan
 #endif

stan/math/prim/scal/meta/partials_return_type.hpp

@@ -25,7 +25,5 @@ namespace stan {
       ::type
       type;
     };
-
 }
 #endif
-

stan/math/prim/scal/meta/partials_type.hpp

@@ -10,4 +10,3 @@ namespace stan {
 
 }
 #endif
-

stan/math/rev/mat/meta/operands_and_partials.hpp

@@ -9,51 +9,31 @@
 namespace stan {
   namespace math {
     namespace detail {
-      template <typename ViewElt, typename Op>
-      class ops_partials_edge_vec_rev
-        : public ops_partials_edge_mat_prim<ViewElt, Op, -1, 1> {
-      public:
-        ops_partials_edge_vec_rev(const Op& ops)
-          : ops_partials_edge_mat_prim<ViewElt, Op, -1, 1>(ops) {}
-        void dump_partials(double* partials) {
-          for (int i = 0; i < this->size(); ++i) {
-            partials[i] = this->partials[i];
-          }
-        }
-        void dump_operands(vari** varis) {
-          for (int i = 0; i < this->size(); ++i) {
-            varis[i] = this->operands[i].vi_;
-          }
-        }
-      };
-
+      // Vectorized Univariate
       // ViewElt = double, Arg = std::vector<var> d_ holds vector<double>
       template <typename ViewElt>
       class ops_partials_edge<ViewElt, std::vector<var> >
-        : public ops_partials_edge_vec_rev<ViewElt, std::vector<var> > {
+        : public ops_partials_edge_mat_prim<ViewElt, std::vector<var>, -1, 1> {
       public:
         explicit ops_partials_edge(const std::vector<var>& ops)
-          : ops_partials_edge_vec_rev<ViewElt, std::vector<var> >(ops) {}
-      };
-
-      // ViewElt = double, Arg = vector_v d_ holds vector<double>
-      // ViewElt = double, Arg = VectorXd d_ holds dummy
-      template <typename ViewElt>
-      class ops_partials_edge<ViewElt, vector_v >
-        : public ops_partials_edge_vec_rev<ViewElt, vector_v > {
-      public:
-        explicit ops_partials_edge(const vector_v& ops)
-          : ops_partials_edge_vec_rev<ViewElt, vector_v>(ops) {}
+          : ops_partials_edge_mat_prim<ViewElt, std::vector<var>, -1, 1>(ops) {}
+        void dump_operands(vari** varis) {
+          for (size_t i = 0; i < this->operands.size(); ++i) {
+            varis[i] = this->operands[i].vi_;
+          }
+        }
       };
-
-      // ViewElt = double, Arg = row_vector_v d_ holds = vector<double>
-      // ViewElt = double, Arg = RowVectorXd, d_ holds dummy
-      template <typename ViewElt>
-      class ops_partials_edge<ViewElt, row_vector_v >
-        : public ops_partials_edge_vec_rev<ViewElt, row_vector_v > {
+      template <typename ViewElt, int R, int C>
+      class ops_partials_edge<ViewElt, Eigen::Matrix<var, R, C> >
+        : public ops_partials_edge_mat_prim<ViewElt, Eigen::Matrix<var, R, C>, R, C> {
       public:
-        explicit ops_partials_edge(const row_vector_v& ops)
-          : ops_partials_edge_vec_rev<ViewElt, row_vector_v>(ops) {}
+        explicit ops_partials_edge(const Eigen::Matrix<var, R, C>& ops)
+          : ops_partials_edge_mat_prim<ViewElt, Eigen::Matrix<var, R, C>, R, C>(ops) {}
+        void dump_operands(vari** varis) {
+          for (int i = 0; i < this->operands.size(); ++i) {
+            varis[i] = this->operands(i).vi_;
+          }
+        }
       };
 
       // MULTIVARIATE
@@ -70,9 +50,9 @@ namespace stan {
           : ops_partials_edge_multivariate_prim<ViewElt, Op>(ops) {}
         void dump_partials(double* partials) {
           int p_i = 0;
-          for (size_t r = 0; r < this->partials.size(); ++r) {
-            for (int c = 0; c < this->partials[r].size(); ++c, ++p_i) {
-              partials[p_i] = this->partials[r](c);
+          for (size_t r = 0; r < this->partials_vec.size(); ++r) {
+            for (int c = 0; c < this->partials_vec[r].size(); ++c, ++p_i) {
+              partials[p_i] = this->partials_vec[r](c);
             }
           }
         }
@@ -102,38 +82,6 @@ namespace stan {
           : ops_partials_edge_multi_rev<ViewElt, vector_v>(ops) {}
       };
 
-      // VIEWS of MATRICES (R = -1, C = -1)
-      // =====================
-      // ViewElt = MatrixXd, Arg = matrix_v, d_ holds MatrixXd
-      // ViewElt = MatrixXd, Arg = MatrixXd d_ holds dummy
-      template <typename ViewElt, int R, int C>
-      class ops_partials_edge<ViewElt, Eigen::Matrix<var, R, C> > {
-      public:
-        typedef Eigen::Matrix<ViewElt, R, C> partials_t;
-        typedef Eigen::Matrix<var, R, C> operands_t;
-        ops_partials_edge(const operands_t& ops)
-          : operands(ops), partials(partials_t::Zero(ops.rows(), ops.cols())) {}
-        void increment_dx_vector(int /*n*/, const partials_t& adj) {
-          partials += adj;
-        }
-        void dump_partials(double* partials) {
-          for (int i = 0; i < this->partials.size(); ++i) {
-            partials[i] = this->partials(i);
-          }
-        }
-        void dump_operands(vari** varis) {
-          for (int i = 0; i < this->operands.size(); ++i) {
-            varis[i] = this->operands(i).vi_;
-          }
-        }
-        int size() {
-          return this->operands.size();
-        }
-      protected:
-        const operands_t& operands;
-        partials_t partials;
-      };
-
       // ViewElt = MatrixXd, Arg = vector<matrix_v> d_ holds vector<MatriXd>
       // ViewElt = MatrixXd, Arg = vector<MatrixXd> d_ holds dummy
       template <typename ViewElt>

stan/math/rev/scal/meta/OperandsAndPartials.hpp

-#ifndef STAN_MATH_REV_SCAL_META_OPERANDSANDPARTIALS_HPP
-#define STAN_MATH_REV_SCAL_META_OPERANDSANDPARTIALS_HPP
-
-#include <stan/math/prim/scal/meta/is_constant_struct.hpp>
-#include <stan/math/prim/scal/meta/length.hpp>
-#include <stan/math/prim/scal/meta/OperandsAndPartials.hpp>
-#include <stan/math/prim/scal/meta/VectorView.hpp>
-#include <stan/math/rev/core.hpp>
-
-namespace stan {
-  namespace math {
-
-    namespace {
-      class partials_vari : public vari {
-      private:
-        const size_t N_;
-        vari** operands_;
-        double* partials_;
-      public:
-        partials_vari(double value,
-                      size_t N,
-                      vari** operands, double* partials)
-          : vari(value),
-            N_(N),
-            operands_(operands),
-            partials_(partials) { }
-        void chain() {
-          for (size_t n = 0; n < N_; ++n)
-            operands_[n]->adj_ += adj_ * partials_[n];
-        }
-      };
-
-      var partials_to_var(double logp, size_t nvaris,
-                          vari** all_varis,
-                          double* all_partials) {
-        return var(new partials_vari(logp, nvaris, all_varis,
-                                     all_partials));
-      }
-
-      template<typename T,
-               bool is_vec = is_vector<T>::value,
-               bool is_const = is_constant_struct<T>::value>
-      struct set_varis {
-        inline size_t set(vari** /*varis*/, const T& /*x*/) {
-          return 0U;
-        }
-      };
-      template<typename T>
-      struct set_varis<T, true, false> {
-        inline size_t set(vari** varis, const T& x) {
-          for (size_t n = 0; n < length(x); n++)
-            varis[n] = x[n].vi_;
-          return length(x);
-        }
-      };
-      template<>
-      struct set_varis<var, false, false> {
-        inline size_t set(vari** varis, const var& x) {
-          varis[0] = x.vi_;
-          return (1);
-        }
-      };
-    }
-
-    /**
-     * This class builds partial derivatives with respect to a set of
-     * operands. There are two reason for the generality of this
-     * class. The first is to handle vector and scalar arguments
-     * without needing to write additional code. The second is to use
-     * this class for writing probability distributions that handle
-     * primitives, reverse mode, and forward mode variables
-     * seamlessly.
-     *
-     * This is the partial template specialization for when the return
-     * type is stan::math::var.
-     *
-     * @tparam T1 First set of operands.
-     * @tparam T2 Second set of operands.
-     * @tparam T3 Third set of operands.
-     * @tparam T4 Fourth set of operands.
-     * @tparam T5 Fifth set of operands.
-     * @tparam T6 Sixth set of operands.
-     * @tparam T_return_type Return type of the expression. This defaults
-     *   to a template metaprogram that calculates the scalar promotion of
-     *   T1 -- T6.
-     */
-    template<typename T1, typename T2, typename T3,
-             typename T4, typename T5, typename T6>
-    struct OperandsAndPartials<T1, T2, T3, T4, T5, T6, stan::math::var> {
-      size_t nvaris;
-      vari** all_varis;
-      double* all_partials;
-
-      VectorView<double,
-                 is_vector<T1>::value,
-                 is_constant_struct<T1>::value> d_x1;
-      VectorView<double,
-                 is_vector<T2>::value,
-                 is_constant_struct<T2>::value> d_x2;
-      VectorView<double,
-                 is_vector<T3>::value,
-                 is_constant_struct<T3>::value> d_x3;
-      VectorView<double,
-                 is_vector<T4>::value,
-                 is_constant_struct<T4>::value> d_x4;
-      VectorView<double,
-                 is_vector<T5>::value,
-                 is_constant_struct<T5>::value> d_x5;
-      VectorView<double,
-                 is_vector<T6>::value,
-                 is_constant_struct<T6>::value> d_x6;
-
-      /**
-       * Constructor.
-       *
-       * @param x1 first set of operands
-       * @param x2 second set of operands
-       * @param x3 third set of operands
-       * @param x4 fourth set of operands
-       * @param x5 fifth set of operands
-       * @param x6 sixth set of operands
-       */
-      OperandsAndPartials(const T1& x1 = 0, const T2& x2 = 0, const T3& x3 = 0,
-                          const T4& x4 = 0, const T5& x5 = 0, const T6& x6 = 0)
-        : nvaris(!is_constant_struct<T1>::value * length(x1) +
-                 !is_constant_struct<T2>::value * length(x2) +
-                 !is_constant_struct<T3>::value * length(x3) +
-                 !is_constant_struct<T4>::value * length(x4) +
-                 !is_constant_struct<T5>::value * length(x5) +
-                 !is_constant_struct<T6>::value * length(x6)),
-          // TODO(carpenter): replace with array allocation fun
-          all_varis(static_cast<vari**>
-                    (vari::operator new
-                     (sizeof(vari*) * nvaris))),
-          all_partials(static_cast<double*>
-                       (vari::operator new
-                        (sizeof(double) * nvaris))),
-          d_x1(all_partials),
-          d_x2(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)),
-          d_x3(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)),
-          d_x4(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)
-               + (!is_constant_struct<T3>::value) * length(x3)),
-          d_x5(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)
-               + (!is_constant_struct<T3>::value) * length(x3)
-               + (!is_constant_struct<T4>::value) * length(x4)),
-          d_x6(all_partials
-               + (!is_constant_struct<T1>::value) * length(x1)
-               + (!is_constant_struct<T2>::value) * length(x2)
-               + (!is_constant_struct<T3>::value) * length(x3)
-               + (!is_constant_struct<T4>::value) * length(x4)
-               + (!is_constant_struct<T5>::value) * length(x5)) {
-        size_t base = 0;
-        if (!is_constant_struct<T1>::value)
-          base += set_varis<T1>().set(&all_varis[base], x1);
-        if (!is_constant_struct<T2>::value)
-          base += set_varis<T2>().set(&all_varis[base], x2);
-        if (!is_constant_struct<T3>::value)
-          base += set_varis<T3>().set(&all_varis[base], x3);
-        if (!is_constant_struct<T4>::value)
-          base += set_varis<T4>().set(&all_varis[base], x4);
-        if (!is_constant_struct<T5>::value)
-          base += set_varis<T5>().set(&all_varis[base], x5);
-        if (!is_constant_struct<T6>::value)
-          set_varis<T6>().set(&all_varis[base], x6);
-        std::fill(all_partials, all_partials+nvaris, 0);
-      }
-
-      /**
-       * Returns a T_return_type with the value specified with
-       * the partial derivatves.
-       *
-       * @param[in] value Value of the variable
-       * @returns a variable with the appropriate value and
-       *   the adjoints set for reverse mode autodiff
-       */
-      stan::math::var value(double value) {
-        return partials_to_var(value, nvaris, all_varis,
-                               all_partials);
-      }
-    };
-
-  }
-}
-#endif

stan/math/rev/scal/meta/operands_and_partials.hpp

@@ -19,83 +19,50 @@ namespace stan {
         ops_partials_edge(const var& op)
           : ops_partials_edge_singular<ViewElt, var>(op) {}
         void dump_partials(ViewElt* partials) {
-          *partials = this->partial;
+          *partials = this->partials[0];
         }
         void dump_operands(vari** varis) {
           *varis = this->operand.vi_;
         }
       };
+    } // end namespace detail
+    template <typename Op1, typename Op2, typename Op3, typename Op4>
+    class operands_and_partials<Op1, Op2, Op3, Op4, var> {
+    public:
+      // these are going to be stack local and get collected
+      detail::ops_partials_edge<double, Op1> edge1_;
+      detail::ops_partials_edge<double, Op2> edge2_;
+      detail::ops_partials_edge<double, Op3> edge3_;
+      detail::ops_partials_edge<double, Op4> edge4_;
+      operands_and_partials(const Op1& o1)
+        : edge1_(o1) { }
+      operands_and_partials(const Op1& o1, const Op2& o2)
+        : edge1_(o1), edge2_(o2) { }
+      operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3)
+        : edge1_(o1), edge2_(o2), edge3_(o3) { }
+      operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3, const Op4& o4)
+        : edge1_(o1), edge2_(o2), edge3_(o3), edge4_(o4) { }
 
-      template <typename Op1, typename Op2, typename Op3, typename Op4>
-      class operands_and_partials<Op1, Op2, Op3, Op4, var> {
-      public:
-        operands_and_partials(const Op1& o1)
-          : edge1_(o1) { }
-        operands_and_partials(const Op1& o1, const Op2& o2)
-          : edge1_(o1), edge2_(o2) { }
-        operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3)
-          : edge1_(o1), edge2_(o2), edge3_(o3) { }
-        operands_and_partials(const Op1& o1, const Op2& o2, const Op3& o3, const Op4& o4)
-          : edge1_(o1), edge2_(o2), edge3_(o3), edge4_(o4) { }
-
-        void increment_dx1(const int n, const double adj) { edge1_.increment_dx(n, adj); }
-        void increment_dx2(const int n, const double adj) { edge2_.increment_dx(n, adj); }
-        void increment_dx3(const int n, const double adj) { edge3_.increment_dx(n, adj); }
-        void increment_dx4(const int n, const double adj) { edge4_.increment_dx(n, adj); }
-
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op1>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx1_vector(const int n, const T& adj) {
-          edge1_.increment_dx_vector(n, adj);
-        }
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op2>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx2_vector(const int n, const T& adj) {
-          edge2_.increment_dx_vector(n, adj);
-        }
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op3>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx3_vector(const int n, const T& adj) {
-          edge3_.increment_dx_vector(n, adj);
-        }
-        template<typename T>
-        typename boost::enable_if_c<is_vector_like<Op4>::value
-                                    && is_vector_like<T>::value, void>::type
-        increment_dx4_vector(const int n, const T& adj) {
-          edge4_.increment_dx_vector(n, adj);
-        }
-
-        // this is what matters in terms of going on autodiff stack
-        var build(double value) {
-          size_t size = edge1_.size() + edge2_.size() + edge3_.size() + edge4_.size();
-          vari** varis = ChainableStack::memalloc_.alloc_array<vari*>(size);
-          int idx = 0;
-          edge1_.dump_operands(&varis[idx]);
-          edge2_.dump_operands(&varis[idx += edge1_.size()]);
-          edge3_.dump_operands(&varis[idx += edge2_.size()]);
-          edge4_.dump_operands(&varis[idx += edge3_.size()]);
+      // this is what matters in terms of going on autodiff stack
+      var build(double value) {
+        size_t size = edge1_.size() + edge2_.size() + edge3_.size() + edge4_.size();
+        vari** varis = ChainableStack::memalloc_.alloc_array<vari*>(size);
+        int idx = 0;
+        edge1_.dump_operands(&varis[idx]);
+        edge2_.dump_operands(&varis[idx += edge1_.size()]);
+        edge3_.dump_operands(&varis[idx += edge2_.size()]);
+        edge4_.dump_operands(&varis[idx += edge3_.size()]);
 
-          double* partials = ChainableStack::memalloc_.alloc_array<double>(size);
-          idx = 0;
-          edge1_.dump_partials(&partials[idx]);
-          edge2_.dump_partials(&partials[idx += edge1_.size()]);
-          edge3_.dump_partials(&partials[idx += edge2_.size()]);
-          edge4_.dump_partials(&partials[idx += edge3_.size()]);
+        double* partials = ChainableStack::memalloc_.alloc_array<double>(size);
+        idx = 0;
+        edge1_.dump_partials(&partials[idx]);
+        edge2_.dump_partials(&partials[idx += edge1_.size()]);
+        edge3_.dump_partials(&partials[idx += edge2_.size()]);
+        edge4_.dump_partials(&partials[idx += edge3_.size()]);
 
-          return var(new precomputed_gradients_vari(value, size, varis, partials));
-        };
-      private:
-        // these are going to be stack local and get collected
-        // TODO(sean): should we pass in ViewElt somehow? shape of it?
-        ops_partials_edge<double, Op1> edge1_;
-        ops_partials_edge<double, Op2> edge2_;
-        ops_partials_edge<double, Op3> edge3_;
-        ops_partials_edge<double, Op4> edge4_;
+        return var(new precomputed_gradients_vari(value, size, varis, partials));
       };
-    } // end namespace detail
+    };
   }
 }
 #endif

test/foo_test.cpp

-#include <stan/math/mix/scal.hpp>
-#include <gtest/gtest.h>
-
-TEST(foo, bar) {
-  double y = 0;
-  double mu = 0;
-  stan::math::fvar<stan::math::var> sigma = 1;
-
-
-  stan::math::fvar<stan::math::var> lp = normal_lpdf(y, mu, sigma);
-}

test/unit/math/fwd/arr/meta/OperandsAndPartials_test.cpp

-#include <stan/math/fwd/arr.hpp>
-#include <gtest/gtest.h>
-
-TEST(AgradPartialsVari, OperandsAndPartialsFvarVec) {
-  using stan::math::OperandsAndPartials;
-  using stan::math::fvar;
-
-  std::vector<fvar<double> > x1;
-  x1.push_back(fvar<double>(2.0,2.0));
-  x1.push_back(fvar<double>(1.0,3.0));
-
-  fvar<double> x2 = 3.0;
-  fvar<double> x3 = 5.0;
-  x2.d_ = -1.0;
-  x3.d_ = 4.0;
-
-  OperandsAndPartials<std::vector<fvar<double> >,fvar<double>,fvar<double> > o(x1, x2, x3);
-  o.d_x1[0] += 17.0; 
-  o.d_x1[1] += 13.0; 
-  o.d_x2[0] += 19.0;  
-  o.d_x2[0] += 19.0;  
-  o.d_x3[0] += 23.0;
-  o.d_x3[0] += 23.0;
-  fvar<double> y = o.value(-1.0);
-
-  EXPECT_FLOAT_EQ(2*17 + 3*13 - 2*19 + 2*4*23,y.d_);
-  EXPECT_FLOAT_EQ(-1,y.val_);
-}
-
-// TEST(AgradPartialsVari, incr_deriv_vec_fvar) {
-//   using stan::VectorView;
-//   using stan::math::incr_deriv;
-//   using stan::is_vector;
-//   using stan::is_constant_struct;
-//   using stan::math::fvar;
-
-//   fvar<double> c(1,1);
-
-//   std::vector<fvar<double> > d;
-//   d.push_back(c);
-//   d.push_back(c);
-  
-//   std::vector<double> d_deriv;
-//   d_deriv.push_back(3);
-//   d_deriv.push_back(4);
-
-//   VectorView<std::vector<double>,
-//              stan::is_vector<std::vector<fvar<double> > >::value,
-//              stan::is_constant_struct<std::vector<fvar<double> > >::value>
-//     d_d(d_deriv);
-
-//   double result3 = incr_deriv<VectorView<std::vector<double>,
-//                                          is_vector<std::vector<fvar<double> > >::value,
-//                                          is_constant_struct<std::vector<fvar<double> > >::value>,
-//                               std::vector<fvar<double> >,double>().incr(d_d,d);
-
-//   EXPECT_FLOAT_EQ(7, result3);
-// }

test/unit/math/fwd/mat/meta/operands_and_partials_test.cpp

@@ -2,7 +2,7 @@
 #include <gtest/gtest.h>
 
 TEST(AgradPartialsVari, OperandsAndPartialsFvarVec) {
-  using stan::math::detail::operands_and_partials;
+  using stan::math::operands_and_partials;
   using stan::math::fvar;
 
   std::vector<fvar<double> > x1;
@@ -18,11 +18,11 @@ TEST(AgradPartialsVari, OperandsAndPartialsFvarVec) {
   dx1 << 17.0, 13.0;
 
   operands_and_partials<std::vector<fvar<double> >,fvar<double>,fvar<double> > o(x1, x2, x3);
-  o.increment_dx1_vector(0, dx1);
-  o.increment_dx2(0, 19.0);
-  o.increment_dx2(0, 19.0);
-  o.increment_dx3(0, 23.0);
-  o.increment_dx3(0, 23.0);
+  o.edge1_.partials_vec[0] += dx1;
+  o.edge2_.partials[0] += 19.0;
+  o.edge2_.partials[0] += 19.0;
+  o.edge3_.partials[0] += 23.0;
+  o.edge3_.partials[0] += 23.0;
   fvar<double> y = o.build(-1.0);
 
   EXPECT_FLOAT_EQ(2*17 + 3*13 - 2*19 + 2*4*23,y.d_);

test/unit/math/fwd/scal/meta/OperandsAndPartials_test.cpp

-#include <stan/math/fwd/scal.hpp>
-#include <gtest/gtest.h>
-
-TEST(AgradPartialsVari, OperandsAndPartialsFvar) {
-  using stan::math::OperandsAndPartials;
-  using stan::math::fvar;
-
-  fvar<double> x1 = 2.0;
-  fvar<double> x2 = 3.0;
-  fvar<double> x3 = 5.0;
-  x1.d_ = 2.0;
-  x2.d_ = -1.0;
-  x3.d_ = 4.0;
-
-  OperandsAndPartials<fvar<double>,fvar<double>,fvar<double> > o(x1, x2, x3);
-  o.d_x1[0] += 17.0; 
-  o.d_x2[0] += 19.0;  
-  o.d_x3[0] += 23.0;
-
-  fvar<double> y = o.value(-1.0);
-  EXPECT_FLOAT_EQ(107,y.d_);
-  EXPECT_FLOAT_EQ(-1,y.val_);
-}
-
-// TEST(AgradPartialsVari, incr_deriv_fvar) {
-//   using stan::VectorView;
-//   using stan::math::incr_deriv;
-//   using stan::is_vector;
-//   using stan::is_constant_struct;
-//   using stan::math::fvar;
-
-//   fvar<double> c;
-//   c.val_ = 1.0;
-//   c.d_ = 1.0;
-
-//   double c_deriv = 2;
-//   VectorView<const double,
-//              stan::is_vector<fvar<double> >::value,
-//              stan::is_constant_struct<fvar<double> >::value>
-//     d_c(c_deriv);
-
-//   double result2 = incr_deriv<VectorView<const double,
-//                                          is_vector<fvar<double> >::value,
-//                                          is_constant_struct<fvar<double> >::value>,
-//                               fvar<double>,double>().incr(d_c,c);
-
-//   EXPECT_FLOAT_EQ(2, result2);
-// }

test/unit/math/fwd/scal/meta/operands_and_partials_test.cpp

@@ -2,7 +2,7 @@
 #include <gtest/gtest.h>
 
 TEST(AgradPartialsVari, OperandsAndPartialsFvar) {
-  using stan::math::detail::operands_and_partials;
+  using stan::math::operands_and_partials;
   using stan::math::fvar;
 
   fvar<double> x1 = 2.0;
@@ -13,9 +13,9 @@ TEST(AgradPartialsVari, OperandsAndPartialsFvar) {
   x3.d_ = 4.0;
 
   operands_and_partials<fvar<double>,fvar<double>,fvar<double> > o(x1, x2, x3);
-  o.increment_dx1(0, 17.0);
-  o.increment_dx2(0, 19.0);
-  o.increment_dx3(0, 23.0);
+  o.edge1_.partials[0] += 17.0;
+  o.edge2_.partials[0] += 19.0;
+  o.edge3_.partials[0] += 23.0;
 
   fvar<double> y = o.build(-1.0);
   EXPECT_FLOAT_EQ(107,y.d_);