Untangles all of our includes!

.gitignore

@@ -19,6 +19,7 @@
 !/test/*.*
 /test/**/*_generated_*_test.cpp
 /test/**/*.xml
+/test/dummy.cpp
 
 # other binary
 *.o

makefile

@@ -109,6 +109,7 @@ endif
 	@echo ''
 	@echo '  To run a single header test, add "-test" to the end of the file name.'
 	@echo '  Example: make stan/math/constants.hpp-test'
+	@echo ''
 	@echo '  - test-math-dependencies : walks through all the header files and indicates'
 	@echo '      when the math dependencies are violated. Dependencies should follow:'
 	@echo '      * rev -> prim'

stan/math/fwd/mat/functor/gradient.hpp

@@ -15,9 +15,8 @@ namespace stan {
      * <p>The functor must implement
      *
      * <code>
-     * stan::math::fvar
-     * operator()(const
-     * Eigen::Matrix<stan::math::var, Eigen::Dynamic, 1>&)
+     * stan::math::fvar<T>
+     * operator()(const Eigen::Matrix<T, Eigen::Dynamic, 1>&)
      * </code>
      *
      * using only operations that are defined for

stan/math/fwd/scal.hpp

@@ -4,6 +4,7 @@
 #include <stan/math/fwd/core.hpp>
 #include <stan/math/fwd/scal/meta/is_fvar.hpp>
 #include <stan/math/fwd/scal/meta/partials_type.hpp>
+#include <stan/math/fwd/scal/meta/OperandsAndPartials.hpp>
 
 #include <stan/math/prim/scal.hpp>
 

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

@@ -0,0 +1,188 @@
+#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 {
+
+    // These are helpers to the OperandsAndPartials specialization for
+    // stan::math::fvar
+    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 stan::math::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 stan::math::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,
+                               typename stan::math::fvar<T_partials_return> > {
+      typedef typename stan::math::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/mix/scal.hpp

@@ -4,10 +4,12 @@
 #include <stan/math/fwd/core.hpp>
 #include <stan/math/fwd/scal/meta/is_fvar.hpp>
 #include <stan/math/fwd/scal/meta/partials_type.hpp>
+#include <stan/math/fwd/scal/meta/OperandsAndPartials.hpp>
 
 #include <stan/math/rev/core.hpp>
 #include <stan/math/rev/scal/meta/is_var.hpp>
 #include <stan/math/rev/scal/meta/partials_type.hpp>
+#include <stan/math/rev/scal/meta/OperandsAndPartials.hpp>
 
 #include <stan/math/prim/scal.hpp>
 #include <stan/math/fwd/scal.hpp>

stan/math/prim/arr.hpp

@@ -4,13 +4,21 @@
 #include <stan/math/prim/arr/meta/container_view.hpp>
 #include <stan/math/prim/arr/meta/get.hpp>
 #include <stan/math/prim/arr/meta/index_type.hpp>
+#include <stan/math/prim/arr/meta/is_constant_struct.hpp>
 #include <stan/math/prim/arr/meta/is_vector.hpp>
 #include <stan/math/prim/arr/meta/length.hpp>
+#include <stan/math/prim/arr/meta/value_type.hpp>
+#include <stan/math/prim/arr/meta/VectorBuilderHelper.hpp>
+#include <stan/math/prim/arr/meta/VectorView.hpp>
+
+#include <stan/math/prim/arr/err/check_ordered.hpp>
 
 #include <stan/math/prim/arr/fun/dist.hpp>
 #include <stan/math/prim/arr/fun/dot.hpp>
 #include <stan/math/prim/arr/fun/dot_self.hpp>
 #include <stan/math/prim/arr/fun/log_sum_exp.hpp>
+#include <stan/math/prim/arr/fun/promote_scalar.hpp>
+#include <stan/math/prim/arr/fun/promote_scalar_type.hpp>
 #include <stan/math/prim/arr/fun/rep_array.hpp>
 #include <stan/math/prim/arr/fun/scaled_add.hpp>
 #include <stan/math/prim/arr/fun/sub.hpp>

stan/math/prim/arr/err/check_ordered.hpp

@@ -0,0 +1,56 @@
+#ifndef STAN_MATH_PRIM_ARR_ERR_CHECK_ORDERED_HPP
+#define STAN_MATH_PRIM_ARR_ERR_CHECK_ORDERED_HPP
+
+#include <stan/math/prim/arr/meta/index_type.hpp>
+#include <stan/math/prim/scal/err/domain_error.hpp>
+#include <stan/math/prim/scal/meta/error_index.hpp>
+#include <sstream>
+#include <string>
+#include <vector>
+
+namespace stan {
+  namespace math {
+
+    /**
+     * Return <code>true</code> if the specified vector is sorted into
+     * strictly increasing order.
+     *
+     * @tparam T_y Type of scalar
+     *
+     * @param function Function name (for error messages)
+     * @param name Variable name (for error messages)
+     * @param y <code>std::vector</code> to test
+     *
+     * @return <code>true</code> if the vector is ordered
+     * @throw <code>std::domain_error</code> if the vector elements are
+     *   not ordered, if there are duplicated
+     *   values, or if any element is <code>NaN</code>.
+     */
+    template <typename T_y>
+    bool check_ordered(const char* function,
+                       const char* name,
+                       const std::vector<T_y>& y) {
+      if (y.size() == 0)
+        return true;
+
+      for (size_t n = 1; n < y.size(); n++) {
+        if (!(y[n] > y[n-1])) {
+          std::ostringstream msg1;
+          msg1 << "is not a valid ordered vector."
+               << " The element at " << stan::error_index::value + n
+               << " is ";
+          std::string msg1_str(msg1.str());
+          std::ostringstream msg2;
+          msg2 << ", but should be greater than the previous element, "
+               << y[n-1];
+          std::string msg2_str(msg2.str());
+          domain_error(function, name, y[n],
+                       msg1_str.c_str(), msg2_str.c_str());
+          return false;
+        }
+      }
+      return true;
+    }
+  }
+}
+#endif

stan/math/prim/arr/fun/promote_scalar.hpp

@@ -0,0 +1,44 @@
+#ifndef STAN_MATH_PRIM_ARR_FUN_PROMOTE_SCALAR_HPP
+#define STAN_MATH_PRIM_ARR_FUN_PROMOTE_SCALAR_HPP
+
+#include <stan/math/prim/scal/fun/promote_scalar.hpp>
+#include <stan/math/prim/scal/fun/promote_scalar_type.hpp>
+#include <stan/math/prim/arr/meta/index_type.hpp>
+#include <vector>
+
+namespace stan {
+  namespace math {
+
+    /**
+     * Struct to hold static function for promoting underlying scalar
+     * types.  This specialization is for standard vector inputs.
+     *
+     * @tparam T return scalar type
+     * @tparam S input type for standard vector elements in static
+     * nested function, which must have an underlying scalar type
+     * assignable to T.
+     */
+    template <typename T, typename S>
+    struct promote_scalar_struct<T, std::vector<S> > {
+      /**
+       * Return the standard vector consisting of the recursive
+       * promotion of the elements of the input standard vector to the
+       * scalar type specified by the return template parameter.
+       *
+       * @param x input standard vector.
+       * @return standard vector with values promoted from input vector.
+       */
+      static std::vector<typename promote_scalar_type<T, S>::type>
+      apply(const std::vector<S>& x) {
+        typedef std::vector<typename promote_scalar_type<T, S>::type> return_t;
+        typedef typename index_type<return_t>::type idx_t;
+        return_t y(x.size());
+        for (idx_t i = 0; i < x.size(); ++i)
+          y[i] = promote_scalar_struct<T, S>::apply(x[i]);
+        return y;
+      }
+    };
+
+  }
+}
+#endif

stan/math/prim/arr/fun/promote_scalar_type.hpp

@@ -0,0 +1,28 @@
+#ifndef STAN_MATH_PRIM_ARR_FUN_PROMOTE_SCALAR_TYPE_HPP
+#define STAN_MATH_PRIM_ARR_FUN_PROMOTE_SCALAR_TYPE_HPP
+
+#include <stan/math/prim/scal/fun/promote_scalar_type.hpp>
+#include <vector>
+
+namespace stan {
+  namespace math {
+
+    /**
+     * Template metaprogram to calculate a type for a container whose
+     * underlying scalar is converted from the second template
+     * parameter type to the first.
+     *
+     * @tparam T result scalar type.
+     * @tparam S input type
+     */
+    template <typename T, typename S>
+    struct promote_scalar_type<T, std::vector<S> > {
+      /**
+       * The promoted type.
+       */
+      typedef std::vector<typename promote_scalar_type<T, S>::type> type;
+    };
+
+  }
+}
+#endif