I’m taking my first baby steps into Stan and it looks great!
One thing I’ve attempted to use recently is vector multiplication.
My original code was with PyStan and it simply threw a “Initialization Error”. I then started simplifying the model and moved to CmdStan to get more verbose exceptions. I landed on some nans in my expressions.
I’ve looked across them and couldn’t figure out where the nans come from, so back to good old print-outs for sanity checking. My code now looks like so:
This compiles but crashes at runtime with:
stan/lib/stan_math/lib/eigen_3.3.9/Eigen/src/Core/Product.h:168: Eigen::ProductImpl<Lhs, Rhs, Option, Eigen::Dense>::Scalar Eigen::ProductImpl<Lhs, Rhs, Option, Eigen::Dense>::coeff(Eigen::Index, Eigen::Index) const [with Lhs = Eigen::Matrix<double, -1, 1>; Rhs = Eigen::Transpose<const Eigen::Matrix<double, -1, 1> >; int Option = 0; Eigen::ProductImpl<Lhs, Rhs, Option, Eigen::Dense>::Scalar = double; Eigen::Index = long int]: Assertion `(Option==LazyProduct) || (this->rows() == 1 && this->cols() == 1)’ failed.
What am I missing in my setup?
I’m trying to perform an outer product here. I think I may have simply assumed that vectors might be reshaped to a Nx1 (or 1xN for row vectors) when needed, but I assume this is not the case?
What’s even weirder, is that if I use the above in an expression (e.g. add 0 to the above), it would not throw an error.
Oh, I didn’t think that was intentional.
It is indeed weird, as outer products are handled by eigen in exactly that way: Eigen: Matrix and vector arithmetic
Always good to have a sanity test. It of course works if I go through with
But something is then still missing, I suppose, between Stan and Eigen. Oh well, at least it clears it up, thanks!
Slightly different version I guess, but here is the relevant eigen source code: Eigen/src/Core/Product.h · master · libeigen / eigen · GitLab
eigen_assert( (Option==LazyProduct) || (this->rows() == 1 && this->cols() == 1) );
Maybe this is a case for @bbbales2?
Edit: where the emphasis should probably be on
Further edit: I think this may explain the behavior in expressions? If it’s part of a lazily evaluated expression, this function never gets called or with option set to the right thing. Not sure about the internals here.
I think this may explain the behavior in expressions?
Yeah looks like it.
Thanks for reporting @idantene. I made an issue: Problem with `print(vec * vec')` · Issue #2435 · stan-dev/math · GitHub
Glad to have helped and be helped at the same time :)