update dirichlet_multinomial

Author: SteveBronder

stan/math/prim/fun/divide.hpp

@@ -43,8 +43,10 @@ inline int divide(int x, int y) {
  */
 template <typename T1, typename T2, require_any_eigen_t<T1, T2>* = nullptr,
           require_all_not_st_var<T1, T2>* = nullptr>
-inline auto divide(const T1& m, const T2& c) {
-  return (as_array_or_scalar(m) / as_array_or_scalar(c)).matrix();
+inline auto divide(T1&& m, T2&& c) {
+  return make_holder([](auto&& m_, auto&& c_) {
+    return (as_array_or_scalar(std::forward<decltype(m_)>(m_)) / as_array_or_scalar(std::forward<decltype(c_)>(c_))).matrix();
+  }, std::forward<T1>(m), std::forward<T2>(c));
 }
 
 }  // namespace math

stan/math/prim/prob/dirichlet_multinomial_lpmf.hpp

@@ -77,7 +77,7 @@ return_type_t<T_prior_size> dirichlet_multinomial_lpmf(
   }
 
   // Need to treat as double otherwise Eigen's array log truncates
-  auto ns_array = as_array_or_scalar(ns).template cast<double>();
+  auto ns_array = as_array_or_scalar(ns).template cast<double>().eval();
   partials_return_t<T_prior_size> a_sum = sum(alpha_val);
   partials_return_t<T_prior_size> lp(0.0);
   if (include_summand<propto>::value) {