From 1cc94324388aa8204ba832a335ae0da384f85d0f Mon Sep 17 00:00:00 2001 From: ariddell Date: Sun, 21 May 2017 14:26:19 -0400 Subject: [PATCH] Compile CVODES functions with each Stan model In 2017, model compilation takes ~40-50 seconds. Compiling CVODES functions with each model adds about 7 seconds. Closes #209. --- pystan/model.py | 41 +++++++++++++++++++++++- pystan/tests/test_cvodes.py | 62 +++++++++++++++++++++++++++++++++++++ 2 files changed, 102 insertions(+), 1 deletion(-) create mode 100644 pystan/tests/test_cvodes.py diff --git a/pystan/model.py b/pystan/model.py index 240762a2..e03f4e4d 100644 --- a/pystan/model.py +++ b/pystan/model.py @@ -269,6 +269,45 @@ def __init__(self, file=None, charset='utf-8', model_name="anon_model", s = template.safe_substitute(model_cppname=self.model_cppname) outfile.write(s) + ## cvodes sources + + # cvodes sources are complied and linked together with the Stan model + # extension module. This is not ideal. In theory, build_clib could be + # used to build a library once and models would be complied and then + # linked with this library. This would save 7 or more seconds from every build. + # But such a strategy is frustrated by the + # lack of ``install_clib`` functionality in Python's distutils. + # + # TODO: numpy provides install_clib functionality, use that. + cvodes_src_path = os.path.join(pystan_dir, 'stan', 'lib', 'stan_math', 'lib', 'cvodes_2.9.0', 'src') + cvodes_sources = [ + os.path.join(cvodes_src_path, 'cvodes', 'cvodea.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodea_io.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_band.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_bandpre.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_bbdpre.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_dense.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_diag.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_direct.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_io.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_sparse.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_spbcgs.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_spgmr.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_spils.c'), + os.path.join(cvodes_src_path, 'cvodes', 'cvodes_sptfqmr.c'), + os.path.join(cvodes_src_path, 'nvec_ser', 'nvector_serial.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_band.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_dense.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_direct.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_iterative.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_math.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_nvector.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_spbcgs.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_spgmr.c'), + os.path.join(cvodes_src_path, 'sundials', 'sundials_sptfqmr.c'), + ] + stan_macros = [ ('BOOST_RESULT_OF_USE_TR1', None), ('BOOST_NO_DECLTYPE', None), @@ -291,7 +330,7 @@ def __init__(self, file=None, charset='utf-8', model_name="anon_model", distutils.log.set_verbosity(verbose) extension = Extension(name=self.module_name, language="c++", - sources=[pyx_file], + sources=[pyx_file] + cvodes_sources, define_macros=stan_macros, include_dirs=include_dirs, extra_compile_args=extra_compile_args) diff --git a/pystan/tests/test_cvodes.py b/pystan/tests/test_cvodes.py new file mode 100644 index 00000000..6d93a833 --- /dev/null +++ b/pystan/tests/test_cvodes.py @@ -0,0 +1,62 @@ +import unittest + +import pystan + + +class TestCVODES(unittest.TestCase): + + def test_cvodes_program(self): + # from integrate_ode_bdf.stan + model_code = """ + functions { + real[] sho(real t, + real[] y, + real[] theta, + real[] x, + int[] x_int) { + real dydt[2]; + dydt[1] = y[2]; + dydt[2] = -y[1] - theta[1] * y[2]; + return dydt; + } + } + data { + int T; + real y0_d[2]; + real t0; + real ts[T]; + real theta_d[1]; + real x[0]; + int x_int[0]; + } + parameters { + real y0_p[2]; + real theta_p[1]; + } + model { + real y_hat[T,2]; + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_d, x, x_int); + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_p, x, x_int); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_d, x, x_int); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_p, x, x_int); + + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_d, x, x_int, 1e-10, 1e-10, 1e8); + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_p, x, x_int, 1e-10, 1e-10, 1e8); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_d, x, x_int, 1e-10, 1e-10, 1e8); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_p, x, x_int, 1e-10, 1e-10, 1e8); + } + generated quantities { + real y_hat[T,2]; + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_d, x, x_int); + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_p, x, x_int); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_d, x, x_int); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_p, x, x_int); + + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_d, x, x_int, 1e-10, 1e-10, 1e8); + y_hat = integrate_ode_bdf(sho, y0_d, t0, ts, theta_p, x, x_int, 1e-10, 1e-10, 1e8); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_d, x, x_int, 1e-10, 1e-10, 1e8); + y_hat = integrate_ode_bdf(sho, y0_p, t0, ts, theta_p, x, x_int, 1e-10, 1e-10, 1e8); + } + """ + model = pystan.StanModel(model_code=model_code, verbose=True) + self.assertIsNotNone(model)