Here's the current dot_product in stan/math/matrix.hpp:
template<int R1,int C1,int R2, int C2>
inline double dot_product(const Eigen::Matrix<double, R1, C1>& v1,
const Eigen::Matrix<double, R2, C2>& v2) {
if (v1.rows() != 1 && v1.cols() != 1)
throw std::invalid_argument("v1 must be a vector");
if (v2.rows() != 1 && v2.cols() != 1)
throw std::invalid_argument("v2 must be a vector");
if (v1.size() != v2.size())
throw std::invalid_argument("v1.size() must equal v2.size()");
return v1.dot(v2);
}
Issue 1
-------
This only compiles if one of R1 and C1 is 1
and one of R2 and C2 is 1, so the dynamic vector checks
are redundant (they'll always pass if the code compiled).
I'm writing the dot_self function to get a bit
more familiar with the templating and wanted to bring
some coding issues to everyone's attention.
Here's the current dot_product in stan/math/matrix.hpp:
template<int R1,int C1,int R2, int C2>
inline double dot_product(const Eigen::Matrix<double, R1, C1>& v1,
const Eigen::Matrix<double, R2, C2>& v2) {
if (v1.rows() != 1 && v1.cols() != 1)
throw std::invalid_argument("v1 must be a vector");
if (v2.rows() != 1 && v2.cols() != 1)
throw std::invalid_argument("v2 must be a vector");
if (v1.size() != v2.size())
throw std::invalid_argument("v1.size() must equal v2.size()");
return v1.dot(v2);
}
Issue 1
-------
This only compiles if one of R1 and C1 is 1
and one of R2 and C2 is 1, so the dynamic vector checks
are redundant (they'll always pass if the code compiled).
The reason is that v1.dot(v2) is only defined if both of
v1 and v2 are Eigen vector template specializations of
Matrix. A type Matrix<T,R,C> is a vector template specialization
if and only if R==1 or C==1.
This means the dynamic checks for row/col size are redundant:
if (v1.rows() != 1 && v1.cols() != 1)
throw std::invalid_argument("v1 must be a vector");
if (v2.rows() != 1 && v2.cols() != 1)
throw std::invalid_argument("v2 must be a vector");
So I removed them.
Issue 2
-------
The dynamic size match test will throw a std::invalid_argument
exception. I'm not sure we want NUTS/HMC to recover form
this. Shouldn't we just stop at this point because it indicates
a programming error? I can't think of any place this would happen
where the situation would be corrected with more sampling.
Issue 3
-------
We're compiling without the -NDEBUG flag, meaning runtime
assertions are turned on.
If we were to get rid of the size check
if (v1.size() != v2.size())
throw std::invalid_argument("v1.size() must equal v2.size()");
then the call
return v1.dot(v2);
will cause assertion failure and crash the program.
I think we want to add the -NDEBUG flag back in for speed
and control the error handling ourselves. For instance,
we don't need to check in v.dot(v) for dot_self()!
I'm writing the dot_self function to get a bit
more familiar with the templating
Interesting. I just got a warning regarding dot_self_vari:
In file included from src/test/agrad/distributions/multivariate/continuous/dirichlet_test.cpp:6:0:
src/stan/agrad/matrix.hpp: In static member function 'static double stan::agrad::dot_self_vari::var_dot_self(const Eigen::Matrix<stan::agrad::var, R, C>&)':
src/stan/agrad/matrix.hpp:498:7: warning: no return statement in function returning non-void [-Wreturn-type]
Indeed, there is no return statement. Could that be causing some of your issues?
I can't believe the other compilers didn't
complain.
- Bob