Solving System of ODE with Eigenvalues/Eigenvectors

Modeling
1

I am modeling a system of ODE that describes the kinetics of a compartmental model. I have successfully used the matrix exponential function in Stan to solve the system for both optimization and sampling. I am now trying to use eigenvalues/eigenvectors to solve the system, extracting the real parts of the complex values returned by the new Stan functions eigenvalues and eigenvectors. Optimization and sampling from the prior work fine with the new method, but I get the following error when I try to sample from the posterior:

Chain 1 iodineModel_2: stan/lib/stan_math/lib/eigen_3.3.9/Eigen/src/Eigenvalues/EigenSolver.h:348: Eigen::EigenSolver<_MatrixType>::EigenvectorsType Eigen::EigenSolver<_MatrixType>::eigenvectors() const [with _MatrixType = Eigen::Matrix<stan::math::var_value, -1, -1, 1, -1, -1>; Eigen::EigenSolver<_MatrixType>::EigenvectorsType = Eigen::Matrix<std::complex<stan::math::var_value >, -1, -1, 1, -1, -1>; typename Eigen::NumTraits::Real = stan::math::var_value]: Assertion `m_eigenvectorsOk && “The eigenvectors have not been computed together with the eigenvalues.”’ failed.
Warning: Chain 1 finished unexpectedly!

There are some rather old non-Stan discussions about this message on the internet that speak to some part of the Eigen library not being implemented yet. Can anyone shed some light on the meaning of this error? Thanks.

2

Since no one appears to know off the cuff what causes this error, I worked up a simple example to illustrate the problem, but unfortunately it still somewhat busy. It attempts to find the estimates of the rate constants and intake of iodine given measurements of iodine in the thyroid and the biokinetic model shown in Figure 2 of Riggs Iodine Model. The MCMC example throws the error mentioned above, but the really interesting thing is that optimization works with a given version of the Stan code but fails if print statements are inserted to aid debugging.

iodineModel_3.stan (3.7 KB)
StanIod-Opt4.R (5.2 KB)
iodineModel_2.stan (3.6 KB)

The attached R and Stan code was run in the following environment:

cmdstan_version()
[1] “2.31.0”
sessionInfo()
R version 4.2.2 Patched (2022-11-10 r83330)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Linux Mint 21.1

Matrix products: default
BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0

locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 LC_MONETARY=en_US.UTF-8
[6] LC_MESSAGES=en_US.UTF-8 LC_PAPER=en_US.UTF-8 LC_NAME=C LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] stats graphics grDevices utils datasets methods base

other attached packages:
[1] cmdstanr_0.5.3

loaded via a namespace (and not attached):
[1] rstudioapi_0.14 knitr_1.41 magrittr_2.0.3 tidyselect_1.2.0 munsell_0.5.0 colorspace_2.0-3
[7] R6_2.5.1 rlang_1.0.6 fansi_1.0.3 dplyr_1.0.10 tools_4.2.2 grid_4.2.2
[13] checkmate_2.1.0 gtable_0.3.1 xfun_0.36 utf8_1.2.2 cli_3.5.0 withr_2.5.0
[19] posterior_1.3.1 abind_1.4-5 tibble_3.1.8 lifecycle_1.0.3 processx_3.8.0 tensorA_0.36.2
[25] farver_2.1.1 ggplot2_3.4.0 ps_1.7.2 vctrs_0.5.1 glue_1.6.2 compiler_4.2.2
[31] pillar_1.8.1 generics_0.1.3 scales_1.2.1 backports_1.4.1 distributional_0.3.1 jsonlite_1.8.4
[37] pkgconfig_2.0.3

2 Likes
3

Maybe worth a look from @WardBrian ?

4

I cannot reproduce the C++ error from above. If I manually run the cmdstan executable with those specified settings it actually samples “fine”, it seems (in the sense that the program terminates normally). Perhaps it is a cmdstanr specific issue?

5

Thanks. So you can’t reproduce the error from the code I provided, right? I will try running it in cmdstan to see if I have the same problems.