Merge pull request #539 from stan-dev/bugfix/0346-copy-to-reference-args

fixes #346, replace copy args with refs

stan/math/fwd/mat/fun/unit_vector_constrain.hpp

@@ -31,25 +31,24 @@ namespace stan {
       for (int k = 0; k < y.size(); ++k)
         unit_vector_y.coeffRef(k).val_ = unit_vector_y_t.coeff(k);
 
-      const T squared_norm = dot_self(y_t);
-      const T norm = sqrt(squared_norm);
-      const T inv_norm = inv(norm);
+      T squared_norm = dot_self(y_t);
+      T norm = sqrt(squared_norm);
+      T inv_norm = inv(norm);
       Matrix<T, Eigen::Dynamic, Eigen::Dynamic> J
         = divide(tcrossprod(y_t),  -norm * squared_norm);
 
       for (int m = 0; m < y.size(); ++m) {
         J.coeffRef(m, m) += inv_norm;
-        for (int k = 0; k < y.size(); ++k) {
+        for (int k = 0; k < y.size(); ++k)
           unit_vector_y.coeffRef(k).d_ = J.coeff(k, m);
-        }
       }
       return unit_vector_y;
     }
 
     template <typename T, int R, int C>
     inline Eigen::Matrix<fvar<T>, R, C>
     unit_vector_constrain(const Eigen::Matrix<fvar<T>, R, C>& y, fvar<T>& lp) {
-      const fvar<T> squared_norm = dot_self(y);
+      fvar<T> squared_norm = dot_self(y);
       lp -= 0.5 * squared_norm;
       return unit_vector_constrain(y);
     }

stan/math/prim/mat/err/check_ldlt_factor.hpp

@@ -11,36 +11,29 @@ namespace stan {
   namespace math {
 
     /**
-     * Check if the argument is a valid
-     * <code>LDLT_factor</code>.
+     * Raise domain error if the specified LDLT factor is invalid.  An
+     * <code>LDLT_factor</code> is invalid if it was constructed from
+     * a matrix that is not positive definite.  The check is that the
+     * <code>success()</code> method returns <code>true</code>.
      *
-     * <code>LDLT_factor</code> can be constructed in an invalid
-     * state, so it must be checked. A invalid <code>LDLT_factor</code>
-     * is constructed from a non positive definite matrix.
-     *
-     * @tparam T Type of scalar
-     * @tparam R Rows of the matrix
-     * @tparam C Columns of the matrix
-     *
-     * @param function Function name (for error messages)
-     * @param name Variable name (for error messages)
-     * @param A <code>LDLT_factor</code> to check for validity.
-     *
-     * @throws <code>std::domain_error</code> the LDLT_factor was
-     *   created improperly (A.success() == false)
+     * @tparam T type of scalar
+     * @tparam R rows of the matrix
+     * @tparam C columns of the matrix
+     * @param[in] function function name for error messages
+     * @param[in] name variable name for error messages
+     * @param[in] A LDLT factor to check for validity
+     * @throws <code>std::domain_error</code> if the LDLT factor is
+     *   invalid.
      */
     template <typename T, int R, int C>
-    inline void check_ldlt_factor(const char* function,
-                                  const char* name,
-                                  LDLT_factor<T, R, C> &A) {
+    inline void check_ldlt_factor(const char* function, const char* name,
+                                  LDLT_factor<T, R, C>& A) {
       if (!A.success()) {
         std::ostringstream msg;
-        msg << "is not positive definite. "
-            << "last conditional variance is ";
+        msg << "is not positive definite.  last conditional variance is ";
         std::string msg_str(msg.str());
-        const T too_small = A.vectorD().tail(1)(0);
-        domain_error(function, name, too_small,
-                     msg_str.c_str(), ".");
+        T too_small = A.vectorD().tail(1)(0);
+        domain_error(function, name, too_small, msg_str.c_str(), ".");
       }
     }
 

stan/math/prim/mat/fun/make_nu.hpp

@@ -8,46 +8,36 @@ namespace stan {
   namespace math {
 
     /**
-     * This function calculates the degrees of freedom for the t
-     * distribution that corresponds to the shape parameter in the
-     * Lewandowski et. al. distribution
+     * Return the degrees of freedom for the t distribution that
+     * corresponds to the shape parameter in the LKJ distribution.
      *
-     * @param eta hyperparameter on (0, inf), eta = 1 <-> correlation
-     * matrix is uniform
-     * @param K number of variables in covariance matrix
+     * @param eta LKJ distribution parameter in (0, inf)
+     * @param K number of variables in correlation matrix
      */
-    template<typename T>
-    const Eigen::Array<T, Eigen::Dynamic, 1>
-    make_nu(const T eta, size_t K) {
+    template <typename T>
+    Eigen::Array<T, Eigen::Dynamic, 1> make_nu(const T& eta, size_t K) {
       using Eigen::Array;
       using Eigen::Dynamic;
       using Eigen::Matrix;
       typedef typename index_type<Matrix<T, Dynamic, 1> >::type size_type;
 
-      Array<T, Dynamic, 1> nu(K * (K - 1) / 2);
-
-      T alpha = eta + (K - 2.0) / 2.0;  // from Lewandowski et. al.
-
       // Best (1978) implies nu = 2 * alpha for the dof in a t
       // distribution that generates a beta variate on (-1, 1)
+      Array<T, Dynamic, 1> nu(K * (K - 1) / 2);
+      T alpha = eta + 0.5 * (K - 2.0);  // from Lewandowski et. al.
       T alpha2 = 2.0 * alpha;
-      for (size_type j = 0; j < (K - 1); j++) {
+      for (size_type j = 0; j < (K - 1); ++j)
         nu(j) = alpha2;
-      }
       size_t counter = K - 1;
-      for (size_type i = 1; i < (K - 1); i++) {
+      for (size_type i = 1; i < (K - 1); ++i) {
         alpha -= 0.5;
         alpha2 = 2.0 * alpha;
-        for (size_type j = i + 1; j < K; j++) {
+        for (size_type j = i + 1; j < K; ++j, ++counter)
           nu(counter) = alpha2;
-          counter++;
-        }
       }
       return nu;
     }
 
   }
-
 }
-
 #endif

stan/math/prim/mat/fun/rank.hpp

@@ -2,53 +2,27 @@
 #define STAN_MATH_PRIM_MAT_FUN_RANK_HPP
 
 #include <stan/math/prim/mat/err/check_range.hpp>
-#include <stan/math/prim/mat/fun/Eigen.hpp>
-#include <vector>
+#include <stan/math/prim/mat/meta/index_type.hpp>
 
 namespace stan {
   namespace math {
 
     /**
      * Return the number of components of v less than v[s].
      *
-     * @tparam T Type of elements.
-     * @param[in] v Input vector.
-     * @param[in] s Position in vector.
-     * @return Number of components of v less than v[s].
+     * @tparam C container type
+     * @param[in] v input vector
+     * @param[in] s position in vector
+     * @return number of components of v less than v[s].
      * @throw std::out_of_range if s is out of range.
      */
-    template <typename T>
-    inline int rank(const std::vector<T> & v, int s) {
-      int size = static_cast<int>(v.size());
-      check_range("rank", "v", size, s);
-      --s;
-      int count(0);
-      T compare(v[s]);
-      for (int i = 0; i < size; ++i)
-        if (v[i] < compare)
-          ++count;
-      return count;
-    }
-
-    /**
-     * Return the number of components of v less than v[s].
-     *
-     * @tparam T Type of elements of the vector.
-     * @param[in] v Input vector.
-     * @param s Index for input vector.
-     * @return Number of components of v less than v[s].
-     * @throw std::out_of_range if s is out of range.
-     */
-    template <typename T, int R, int C>
-    inline int rank(const Eigen::Matrix<T, R, C> & v, int s) {
-      int size = v.size();
-      check_range("rank", "v", size, s);
-      --s;
-      const T * vv = v.data();
-      int count(0);
-      T compare(vv[s]);
-      for (int i = 0; i < size; ++i)
-        if (vv[i] < compare)
+    template <typename C>
+    inline int rank(const C& v, int s) {
+      check_range("rank", "v", v.size(), s);
+      --s;  // adjust for indexing by one
+      int count = 0;
+      for (typename index_type<C>::type i = 0; i < v.size(); ++i)
+        if (v[i] < v[s])
           ++count;
       return count;
     }

stan/math/prim/mat/fun/unit_vector_constrain.hpp

@@ -14,19 +14,22 @@ namespace stan {
 
     /**
      * Return the unit length vector corresponding to the free vector y.
-     * See https://en.wikipedia.org/wiki/N-sphere#Generating_random_points
      *
+     * See <a
+     * href="https://en.wikipedia.org/wiki/N-sphere#Generating_random_points">the
+     * Wikipedia page on generating random points on an N-sphere</a>.
+     *
+     * @tparam T Scalar type.
      * @param y vector of K unrestricted variables
      * @return Unit length vector of dimension K
-     * @tparam T Scalar type.
-     **/
+     */
     template <typename T, int R, int C>
     Eigen::Matrix<T, R, C>
     unit_vector_constrain(const Eigen::Matrix<T, R, C>& y) {
       using std::sqrt;
       check_vector("unit_vector_constrain", "y", y);
       check_nonzero_size("unit_vector_constrain", "y", y);
-      const T SN = dot_self(y);
+      T SN = dot_self(y);
       check_positive_finite("unit_vector_constrain", "norm", SN);
       return y / sqrt(SN);
     }
@@ -39,21 +42,19 @@ namespace stan {
      * @return Unit length vector of dimension K
      * @param lp Log probability reference to increment.
      * @tparam T Scalar type.
-     **/
+     */
     template <typename T, int R, int C>
     Eigen::Matrix<T, R, C>
     unit_vector_constrain(const Eigen::Matrix<T, R, C>& y, T& lp) {
       using std::sqrt;
       check_vector("unit_vector_constrain", "y", y);
       check_nonzero_size("unit_vector_constrain", "y", y);
-      const T SN = dot_self(y);
+      T SN = dot_self(y);
       check_positive_finite("unit_vector_constrain", "norm", SN);
       lp -= 0.5 * SN;
       return y / sqrt(SN);
     }
 
   }
-
 }
-
 #endif

stan/math/prim/scal/fun/as_bool.hpp

@@ -5,13 +5,15 @@ namespace stan {
   namespace math {
 
     /**
-     * Return 1 if the argument is unequal to zero and 0 otherwise.
+     * Return true if the argument is not equal to zero (in the
+     * <code>!=</code> operator sense) and false otherwise.
      *
-     * @param x Value.
-     * @return 1 if argument is equal to zero (or NaN) and 0 otherwise.
+     * @tparam T type of scalar
+     * @param x value
+     * @return true if value is not equal to zero
      */
     template <typename T>
-    inline bool as_bool(const T x) {
+    inline bool as_bool(const T& x) {
       return x != 0;
     }
 

stan/math/prim/scal/fun/binary_log_loss.hpp

@@ -1,7 +1,7 @@
 #ifndef STAN_MATH_PRIM_SCAL_FUN_BINARY_LOG_LOSS_HPP
 #define STAN_MATH_PRIM_SCAL_FUN_BINARY_LOG_LOSS_HPP
 
-#include <boost/math/tools/promotion.hpp>
+#include <stan/math/prim/scal/fun/log1m.hpp>
 
 namespace stan {
   namespace math {
@@ -17,15 +17,15 @@ namespace stan {
      *
      * \f$\mbox{logloss}(0, \hat{y}) = -\log (1 - \hat{y}) \f$.
      *
-     * @param y Reference value in { 0 , 1 }.
-     * @param y_hat Response value in [0, 1].
-     * @return Log loss for response given reference value.
+     * @tparam T value type
+     * @param[in] y reference value, either 0 or 1
+     * @param[in] y_hat response value in [0, 1]
+     * @return Log loss for response given reference value
      */
     template <typename T>
-    inline typename boost::math::tools::promote_args<T>::type
-    binary_log_loss(int y, const T y_hat) {
+    inline T binary_log_loss(int y, const T& y_hat) {
       using std::log;
-      return -log(y ? y_hat : (1.0 - y_hat));
+      return y ? -log(y_hat) : -log1m(y_hat);
     }
 
   }

stan/math/prim/scal/fun/corr_constrain.hpp

@@ -2,6 +2,7 @@
 #define STAN_MATH_PRIM_SCAL_FUN_CORR_CONSTRAIN_HPP
 
 #include <stan/math/prim/scal/fun/log1m.hpp>
+#include <stan/math/prim/scal/fun/square.hpp>
 #include <cmath>
 
 namespace stan {
@@ -15,13 +16,13 @@ namespace stan {
      *
      * <p>\f$f(x) = \tanh x = \frac{\exp(2x) - 1}{\exp(2x) + 1}\f$.
      *
-     * @param x Scalar input.
-     * @return Result of transforming the input to fall between -1 and 1.
-     * @tparam T Type of scalar.
+     * @tparam T type of value
+     * @param[in] x value
+     * @return tanh transform of value
      */
     template <typename T>
     inline
-    T corr_constrain(const T x) {
+    T corr_constrain(const T& x) {
       return tanh(x);
     }
 
@@ -36,17 +37,17 @@ namespace stan {
      * <p>\f$\log | \frac{d}{dx} \tanh x  | = \log (1 - \tanh^2 x)\f$.
      *
      * @tparam T Type of scalar.
+     * @param[in] x value
+     * @param[in,out] lp log density accumulator
      */
     template <typename T>
     inline
-    T corr_constrain(const T x, T& lp) {
+    T corr_constrain(const T& x, T& lp) {
       T tanh_x = tanh(x);
-      lp += log1m(tanh_x * tanh_x);
+      lp += log1m(square(tanh_x));
       return tanh_x;
     }
 
   }
-
 }
-
 #endif

stan/math/prim/scal/fun/corr_free.hpp

@@ -19,15 +19,14 @@ namespace stan {
      *          = \mbox{atanh}\, y
      *          = \frac{1}{2} \log \frac{y + 1}{y - 1}\f$.
      *
-     * @param y Correlation scalar input.
-     * @return Free scalar that transforms to the specified input.
-     * @tparam T Type of scalar.
+     * @tparam T Type of correlation
+     * @param[in] y correlation
+     * @return free scalar that transforms to the specified input
      */
     template <typename T>
     inline
-    T corr_free(const T y) {
-      check_bounded("lub_free",
-                    "Correlation variable", y, -1.0, 1.0);
+    T corr_free(const T& y) {
+      check_bounded("lub_free", "Correlation variable", y, -1.0, 1.0);
       return atanh(y);
     }