I’d like to solve a Lyapanov equation where AX + XA’= Q within my stan model.
I have a c++ function using Eigen that can do this, and I have tried to include it in the stan model, but I just get the seemingly standard " Compilation ERROR, function(s)/method(s) not created!" message. My c++ is not great, I know the code works within R / Rcppeigen, but I’m also not at all sure whether Stan will allow this sort of more elaborate / library dependent c++ , and if so, which if any includes are needed.
Ok, the non created function seems to be due to Eigen::internal::matrix_function_solve_triangular_sylvester
Does the fact that it’s an unsupported function in Eigen mean that it can’t be found? Including: #include <unsupported/Eigen/MatrixFunctions>
didn’t help.
edited:
I couldn’t work it out so just copied / adapted that function. I’m now stuck on the autograd, not sure how to pass the vars through… code below works but all the commented out lines are needed and things break.
Any hand holding would be appreciated!
template <typename T1>
Eigen::Matrix<T1, Eigen::Dynamic, Eigen::Dynamic> sylcpp(const Eigen::Matrix<T1, Eigen::Dynamic, Eigen::Dynamic> & A,
const Eigen::Matrix<T1, Eigen::Dynamic, Eigen::Dynamic> & C){
using matrix_t = Eigen::Matrix<T1, Eigen::Dynamic, Eigen::Dynamic>;
matrix_t B = A.transpose();
// Eigen::ComplexSchur<matrix_t> SchurA(A);
// Eigen::MatrixXcd R = SchurA.matrixT();
// Eigen::MatrixXcd U = SchurA.matrixU();
// Eigen::ComplexSchur<Eigen::MatrixXd> SchurB(B);
// Eigen::MatrixXcd S = SchurB.matrixT();
// Eigen::MatrixXcd V = SchurB.matrixU();
// Eigen::MatrixXcd F = (U.adjoint() * C) * V;
// Eigen::MatrixXcd Y = matrix_function_solve_triangular_sylvester(R, S, F);
// Eigen::MatrixXcd X = (U * Y) * V.adjoint();
// return X.real();
return(A*C+C*A);
}
template <typename T0__, typename T1__>
Eigen::Matrix<typename boost::math::tools::promote_args<T0__, T1__>::type, Eigen::Dynamic, Eigen::Dynamic>
syl(const Eigen::Matrix<T0__, Eigen::Dynamic, Eigen::Dynamic>& A,
const Eigen::Matrix<T1__, Eigen::Dynamic, Eigen::Dynamic>& C, std::ostream* pstream__){
Eigen::Matrix<T0__, Eigen::Dynamic,Eigen::Dynamic> X = sylcpp(A,C);
return(X);
}
Thanks. The blogpost code didn’t work for me, maybe I misunderstood something, but the code above is working fine up to simple functions like in the blog, but I can’t manage to pass the var through all the eigen functions. Possibly a complex number issue in there as well and I will have to find a way to do it with reals.