Call GNU Scientific Library from stan

Modeling
21

Hi James,

Thanks alot for help me! I run the model as follows:
mc2 <-
’
functions { real integral_GK(real ll, real ul, real nu); }
model {} // use the integral_GK() function somehow
’

model_tmp2<-stan_model(model_code = mc2, model_name = “external”, allow_undefined = TRUE,
includes = paste0(’\n#include “’,
file.path(getwd(), ‘exp_GK.hpp’), '”\n’))

and the exp_GK.hpp file is:

template <typename T0__, typename T1__,typename T2__>
typename boost::math::tools::promote_args<T0__, T1__,T2__>::type
integral_GK(const T0__& ll, const T1__& ul,const T2__& nu, std::ostream* pstream__) {
boost::math::quadrature::gauss_kronrod<double, 15> integrator;
auto f1 = [&](double x) {
return std::exp(x);
};
double error;
double ll_new=ll * nu;
double ul_new=ul * nu;
double Q = integrator.integrate(f1, ll_new, ul_new, 0, 0, &error)/nu;
return Q;
}

and the error is:
In file included from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/BH/include/boost/config.hpp:39:0,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/BH/include/boost/math/tools/config.hpp:13,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/core/var.hpp:7,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/core/gevv_vvv_vari.hpp:5,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/core.hpp:12,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/mat.hpp:4,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math.hpp:4,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/src/stan/model/model_header.hpp:4,
from filec2060e74890.cpp:8:
C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/BH/include/boost/config/compiler/gcc.hpp:186:0: warning: “BOOST_NO_CXX11_RVALUE_REFERENCES” redefined

define BOOST_NO_CXX11_RVALUE_REFERENCES

^
:0:0: note: this is the location of the previous definition
In file included from filec2060e74890.cpp:41:0:
C:/Users/Nhat Le Thanh Hoang/Documents/exp_GK.hpp: In function ‘typename boost::math::tools::promote_args<RT1, RT2, A>::type modelc204f092276_external_namespace::integral_GK(const T0__&, const T1__&, const T2__&, std::ostream*)’:
C:/Users/Nhat Le Thanh Hoang/Documents/exp_GK.hpp:4:16: error: ‘boost::math::quadrature’ has not been declared
boost::math::quadrature::gauss_kronrod<double, 15> integrator;
^
C:/Users/Nhat Le Thanh Hoang/Documents/exp_GK.hpp:4:42: error: expected primary-expression before ‘double’
boost::math::quadrature::gauss_kronrod<double, 15> integrator;
^
C:/Users/Nhat Le Thanh Hoang/Documents/exp_GK.hpp:5:8: error: ‘f1’ does not name a type
auto f1 = [&](double x) {
^
C:/Users/Nhat Le Thanh Hoang/Documents/exp_GK.hpp:11:14: error: ‘integrator’ was not declared in this scope
double Q = integrator.integrate(f1, ll_new, ul_new, 0, 0, &error)/nu;
^
C:/Users/Nhat Le Thanh Hoang/Documents/exp_GK.hpp:11:35: error: ‘f1’ was not declared in this scope
double Q = integrator.integrate(f1, ll_new, ul_new, 0, 0, &error)/nu;
^
make: *** [filec2060e74890.o] Error 1

ERROR(s) during compilation: source code errors or compiler configuration errors!

Program source:
1:
2: // includes from the plugin
3:
4:
5: // user includes
6: #define STAN__SERVICES__COMMAND_HPP// Code generated by Stan version 2.17.0
7:
8: #include <stan/model/model_header.hpp>
9:
10: namespace modelc204f092276_external_namespace {
11:
12: using std::istream;
13: using std::string;
14: using std::stringstream;
15: using std::vector;
16: using stan::io::dump;
17: using stan::math::lgamma;
18: using stan::model::prob_grad;
19: using namespace stan::math;
20:
21: typedef Eigen::Matrix<double,Eigen::Dynamic,1> vector_d;
22: typedef Eigen::Matrix<double,1,Eigen::Dynamic> row_vector_d;
23: typedef Eigen::Matrix<double,Eigen::Dynamic,Eigen::Dynamic> matrix_d;
24:
25: static int current_statement_begin__;
26:
27: stan::io::program_reader prog_reader__() {
28: stan::io::program_reader reader;
29: reader.add_event(0, 0, “start”, “modelc204f092276_external”);
30: reader.add_event(3, 3, “end”, “modelc204f092276_external”);
31: return reader;
32: }
33:
34: template <typename T0__, typename T1__, typename T2__>
35: typename boost::math::tools::promote_args<T0__, T1__, T2__>::type
36: integral_GK(const T0__& ll,
37: const T1__& ul,
38: const T2__& nu, std::ostream* pstream__);
39:
40:
41: #include “C:/Users/Nhat Le Thanh Hoang/Documents/exp_GK.hpp”
42: class modelc204f092276_external : public prob_grad {
43: private:
44: public:
45: modelc204f092276_external(stan::io::var_context& context__,
46: std::ostream* pstream__ = 0)
47: : prob_grad(0) {
48: ctor_body(context__, 0, pstream__);
49: }
50:
51: modelc204f092276_external(stan::io::var_context& context__,
52: unsigned int random_seed__,
53: std::ostream* pstream__ = 0)
54: : prob_grad(0) {
55: ctor_body(context__, random_seed__, pstream__);
56: }
57:
58: void ctor_body(stan::io::var_context& context__,
59: unsigned int random_seed__,
60: std::ostream* pstream__) {
61: boost::ecuyer1988 base_rng__ =
62: stan::services::util::create_rng(random_seed__, 0);
63: (void) base_rng__; // suppress unused var warning
64:
65: current_statement_begin__ = -1;
66:
67: static const char* function__ = “modelc204f092276_external_namespace::modelc204f092276_external”;
68: (void) function__; // dummy to suppress unused var warning
69: size_t pos__;
70: (void) pos__; // dummy to suppress unused var warning
71: std::vector vals_i__;
72: std::vector vals_r__;
73: double DUMMY_VAR__(std::numeric_limits::quiet_NaN());
74: (void) DUMMY_VAR__; // suppress unused var warning
75:
76: // initialize member variables
77: try {
78:
79: // validate, data variables
80: // initialize data variables
81:
82:
83: // validate transformed data
84:
85: // validate, set parameter ranges
86: num_params_r__ = 0U;
87: param_ranges_i__.clear();
88: } catch (const std::exception& e) {
89: stan::lang::rethrow_located(e, current_statement_begin__, prog_reader__());
90: // Next line prevents compiler griping about no return
91: throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ");
92: }
93: }
94:
95: ~modelc204f092276_external() { }
96:
97:
98: void transform_inits(const stan::io::var_context& context__,
99: std::vector& params_i__,
100: std::vector& params_r__,
101: std::ostream pstream__) const {
102: stan::io::writer writer__(params_r__,params_i__);
103: size_t pos__;
104: (void) pos__; // dummy call to supress warning
105: std::vector vals_r__;
106: std::vector vals_i__;
107:
108: params_r__ = writer__.data_r();
109: params_i__ = writer__.data_i();
110: }
111:
112: void transform_inits(const stan::io::var_context& context,
113: Eigen::Matrix<double,Eigen::Dynamic,1>& params_r,
114: std::ostream pstream__) const {
115: std::vector params_r_vec;
116: std::vector params_i_vec;
117: transform_inits(context, params_i_vec, params_r_vec, pstream__);
118: params_r.resize(params_r_vec.size());
119: for (int i = 0; i < params_r.size(); ++i)
120: params_r(i) = params_r_vec[i];
121: }
122:
123:
124: template <bool propto__, bool jacobian__, typename T__>
125: T__ log_prob(vector<T__>& params_r__,
126: vector& params_i__,
127: std::ostream pstream__ = 0) const {
128:
129: T__ DUMMY_VAR__(std::numeric_limits::quiet_NaN());
130: (void) DUMMY_VAR__; // suppress unused var warning
131:
132: T__ lp__(0.0);
133: stan::math::accumulator<T__> lp_accum__;
134:
135: try {
136: // model parameters
137: stan::io::reader<T__> in__(params_r__,params_i__);
138:
139:
140: // transformed parameters
141:
142:
143:
144: // validate transformed parameters
145:
146: const char* function__ = “validate transformed params”;
147: (void) function__; // dummy to suppress unused var warning
148:
149: // model body
150:
151:
152: } catch (const std::exception& e) {
153: stan::lang::rethrow_located(e, current_statement_begin__, prog_reader__());
154: // Next line prevents compiler griping about no return
155: throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ");
156: }
157:
158: lp_accum__.add(lp__);
159: return lp_accum__.sum();
160:
161: } // log_prob()
162:
163: template <bool propto, bool jacobian, typename T_>
164: T_ log_prob(Eigen::Matrix<T_,Eigen::Dynamic,1>& params_r,
165: std::ostream pstream = 0) const {
166: std::vector<T_> vec_params_r;
167: vec_params_r.reserve(params_r.size());
168: for (int i = 0; i < params_r.size(); ++i)
169: vec_params_r.push_back(params_r(i));
170: std::vector vec_params_i;
171: return log_prob<propto,jacobian,T_>(vec_params_r, vec_params_i, pstream);
172: }
173:
174:
175: void get_param_names(std::vectorstd::string& names__) const {
176: names__.resize(0);
177: }
178:
179:
180: void get_dims(std::vector<std::vector<size_t> >& dimss__) const {
181: dimss__.resize(0);
182: std::vector<size_t> dims__;
183: }
184:
185: template
186: void write_array(RNG& base_rng__,
187: std::vector& params_r__,
188: std::vector& params_i__,
189: std::vector& vars__,
190: bool include_tparams__ = true,
191: bool include_gqs__ = true,
192: std::ostream pstream__ = 0) const {
193: vars__.resize(0);
194: stan::io::reader in__(params_r__,params_i__);
195: static const char function__ = “modelc204f092276_external_namespace::write_array”;
196: (void) function__; // dummy to suppress unused var warning
197: // read-transform, write parameters
198:
199: if (!include_tparams__) return;
200: // declare and define transformed parameters
201: double lp__ = 0.0;
202: (void) lp__; // dummy to suppress unused var warning
203: stan::math::accumulator lp_accum__;
204:
205: double DUMMY_VAR__(std::numeric_limits::quiet_NaN());
206: (void) DUMMY_VAR__; // suppress unused var warning
207:
208: try {
209:
210:
211:
212: // validate transformed parameters
213:
214: // write transformed parameters
215:
216: if (!include_gqs__) return;
217: // declare and define generated quantities
218:
219:
220:
221: // validate generated quantities
222:
223: // write generated quantities
224: } catch (const std::exception& e) {
225: stan::lang::rethrow_located(e, current_statement_begin__, prog_reader__());
226: // Next line prevents compiler griping about no return
227: throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ");
228: }
229: }
230:
231: template
232: void write_array(RNG& base_rng,
233: Eigen::Matrix<double,Eigen::Dynamic,1>& params_r,
234: Eigen::Matrix<double,Eigen::Dynamic,1>& vars,
235: bool include_tparams = true,
236: bool include_gqs = true,
237: std::ostream pstream = 0) const {
238: std::vector params_r_vec(params_r.size());
239: for (int i = 0; i < params_r.size(); ++i)
240: params_r_vec[i] = params_r(i);
241: std::vector vars_vec;
242: std::vector params_i_vec;
243: write_array(base_rng,params_r_vec,params_i_vec,vars_vec,include_tparams,include_gqs,pstream);
244: vars.resize(vars_vec.size());
245: for (int i = 0; i < vars.size(); ++i)
246: vars(i) = vars_vec[i];
247: }
248:
249: static std::string model_name() {
250: return “modelc204f092276_external”;
251: }
252:
253:
254: void constrained_param_names(std::vectorstd::string& param_names__,
255: bool include_tparams__ = true,
256: bool include_gqs__ = true) const {
257: std::stringstream param_name_stream__;
258:
259: if (!include_gqs__ && !include_tparams__) return;
260:
261: if (!include_gqs__) return;
262: }
263:
264:
265: void unconstrained_param_names(std::vectorstd::string& param_names__,
266: bool include_tparams__ = true,
267: bool include_gqs__ = true) const {
268: std::stringstream param_name_stream__;
269:
270: if (!include_gqs__ && !include_tparams__) return;
271:
272: if (!include_gqs__) return;
273: }
274:
275: }; // model
276:
277: }
278:
279: typedef modelc204f092276_external_namespace::modelc204f092276_external stan_model;
280:
281: #include <rstan/rstaninc.hpp>
282: /
283: * Define Rcpp Module to expose stan_fit’s functions to R.
284: */
285: RCPP_MODULE(stan_fit4modelc204f092276_external_mod){
286: Rcpp::class_<rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external,
287: boost::random::ecuyer1988> >(“stan_fit4modelc204f092276_external”)
288: // .constructorRcpp::List()
289: .constructor<SEXP, SEXP, SEXP>()
290: // .constructor<SEXP, SEXP>()
291: .method(“call_sampler”,
292: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::call_sampler)
293: .method(“param_names”,
294: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::param_names)
295: .method(“param_names_oi”,
296: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::param_names_oi)
297: .method(“param_fnames_oi”,
298: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::param_fnames_oi)
299: .method(“param_dims”,
300: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::param_dims)
301: .method(“param_dims_oi”,
302: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::param_dims_oi)
303: .method(“update_param_oi”,
304: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::update_param_oi)
305: .method(“param_oi_tidx”,
306: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::param_oi_tidx)
307: .method(“grad_log_prob”,
308: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::grad_log_prob)
309: .method(“log_prob”,
310: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::log_prob)
311: .method(“unconstrain_pars”,
312: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::unconstrain_pars)
313: .method(“constrain_pars”,
314: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::constrain_pars)
315: .method(“num_pars_unconstrained”,
316: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::num_pars_unconstrained)
317: .method(“unconstrained_param_names”,
318: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::unconstrained_param_names)
319: .method(“constrained_param_names”,
320: &rstan::stan_fit<modelc204f092276_external_namespace::modelc204f092276_external, boost::random::ecuyer1988>::constrained_param_names)
321: ;
322: }
323:
324: // declarations
325: extern “C” {
326: SEXP filec2060e74890( ) ;
327: }
328:
329: // definition
330:
331: SEXP filec2060e74890( ){
332: return Rcpp::wrap(“external”);
333: }
334:
335:
Error in compileCode(f, code, language = language, verbose = verbose) :
Compilation ERROR, function(s)/method(s) not created! In file included from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/BH/include/boost/config.hpp:39:0,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/BH/include/boost/math/tools/config.hpp:13,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/core/var.hpp:7,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/core/gevv_vvv_vari.hpp:5,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/core.hpp:12,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math/rev/mat.hpp:4,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-library/3.5/StanHeaders/include/stan/math.hpp:4,
from C:/Users/Nhat Le Thanh Hoang/Documents/R/win-
In addition: Warning message:
In system(cmd, intern = !verbose) :
running command ‘C:/PROGRA~1/R/R-35~1.0/bin/x64/R CMD SHLIB filec2060e74890.cpp 2> filec2060e74890.cpp.err.txt’ had status 1

Nhat

22

Offhand, I think this error is due to not having #include <boost/math/quadrature/gauss_kronrod.hpp> at the top of the file exp_GK.hpp, since the compiler doesn’t recognize boost::math::quadrature.

This error is almost certainly due to RStan not using a compiler that accepts C++11 by default. See here: auto - C++ Forum

I’m not entirely sure how to solve this, but judging from First steps in using C++11 with Rcpp, it looks like there are two possible solutions, either adding the comment // [[Rcpp::plugins(cpp11)]] to the top of of the file exp_GK.hpp or running the command Sys.setenv("PKG_CXXFLAGS"="-std=c++11") before invoking stan_model. I’d try the comment first to see what happens.

23

Hi James,

Adding the comment // [[Rcpp::plugins(cpp11)]] to the top of the file make the code works however I’m still unable to include the #include <boost/math/quadrature/gauss_kronrod.hpp> on the top of the file.
Accordding to the instruction https://cran.r-project.org/web/packages/rstan/vignettes/external.html in the example of besselK function. this function in the library boost::math whose header are pulled in the stan math library. The fact that I cannot call the class gauss_kronrod of the library boost/math/quadrature is due to the header are not pulled in the stan math library.

How can I fix this?

Thanks alot!
Nhat

24

Hi all,
I think I have to stop trying to use the class gauss_kronrod in Boost C++ library. It seems out of my hand. I check all the special function in stan and I found that there is a function expint helps us to compute the En integral https://www.boost.org/doc/libs/1_67_0/boost/math/special_functions/expint.hpp
I run the code exp_int.hpp code as follows
template
typename boost::math::tools::promote_args<T0__>::result_type
integrate(unsigned n, const T0__& z, std::ostream* pstream__) {
return boost::math::expint(n, z);
}

and the model
mc <-
’
functions {
real integrate(int n, real z);
}
parameters {real x;}
model {x ~ normal(0,1);}
generated quantities{
real y;

y= x+integrate(0,-1.0);

} // use the exp_int() function somehow
’

model_tmp2<-stan_model(model_code = mc, model_name = “external”, allow_undefined = TRUE,
includes = paste0(’\n#include “’,
file.path(getwd(), ‘exp_int.hpp’), '”\n’))

I have the error as follows:

In file included from C:/Users/nhatlth/Documents/R/win-library/3.3/BH/include/boost/config.hpp:39:0,
from C:/Users/nhatlth/Documents/R/win-library/3.3/BH/include/boost/math/tools/config.hpp:13,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/core/var.hpp:7,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/core/gevv_vvv_vari.hpp:5,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/core.hpp:12,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/mat.hpp:4,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math.hpp:4,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/src/stan/model/model_header.hpp:4,
from file29dc94a1a87.cpp:8:
C:/Users/nhatlth/Documents/R/win-library/3.3/BH/include/boost/config/compiler/gcc.hpp:186:0: warning: “BOOST_NO_CXX11_RVALUE_REFERENCES” redefined

define BOOST_NO_CXX11_RVALUE_REFERENCES

^
:0:0: note: this is the location of the previous definition
file29dc94a1a87.o:file29dc94a1a87.cpp:(.text+0xade5): undefined reference to boost::math::tools::promote_args<double, float, float, float, float, float>::type model29dc489d44f4_external_namespace::integrate<double>(int const&, double const&, std::ostream*)' file29dc94a1a87.o:file29dc94a1a87.cpp:(.text$_ZNK36model29dc489d44f4_external_namespace26model29dc489d44f4_external11write_arrayIN5boost6random23additive_combine_engineINS3_26linear_congruential_engineIjLj40014ELj0ELj2147483563EEENS5_IjLj40692ELj0ELj2147483399EEEEEEEvRT_RSt6vectorIdSaIdEERSB_IiSaIiEESE_bbPSo[_ZNK36model29dc489d44f4_external_namespace26model29dc489d44f4_external11write_arrayIN5boost6random23additive_combine_engineINS3_26linear_congruential_engineIjLj40014ELj0ELj2147483563EEENS5_IjLj40692ELj0ELj2147483399EEEEEEEvRT_RSt6vectorIdSaIdEERSB_IiSaIiEESE_bbPSo]+0xd9): undefined reference toboost::math::tools::promote_args<double, float, float, float, float, float>::type model29dc489d44f4_external_namespace::integrate(int const&, double const&, std::ostream*)’
collect2.exe: error: ld returned 1 exit status

ERROR(s) during compilation: source code errors or compiler configuration errors!

Program source:
1:
2: // includes from the plugin
3:
4:
5: // user includes
6: #define STAN__SERVICES__COMMAND_HPP// Code generated by Stan version 2.17.0
7:
8: #include <stan/model/model_header.hpp>
9:
10: namespace model29dc489d44f4_external_namespace {
11:
12: using std::istream;
13: using std::string;
14: using std::stringstream;
15: using std::vector;
16: using stan::io::dump;
17: using stan::math::lgamma;
18: using stan::model::prob_grad;
19: using namespace stan::math;
20:
21: typedef Eigen::Matrix<double,Eigen::Dynamic,1> vector_d;
22: typedef Eigen::Matrix<double,1,Eigen::Dynamic> row_vector_d;
23: typedef Eigen::Matrix<double,Eigen::Dynamic,Eigen::Dynamic> matrix_d;
24:
25: static int current_statement_begin__;
26:
27: stan::io::program_reader prog_reader__() {
28: stan::io::program_reader reader;
29: reader.add_event(0, 0, “start”, “model29dc489d44f4_external”);
30: reader.add_event(11, 11, “end”, “model29dc489d44f4_external”);
31: return reader;
32: }
33:
34: template
35: typename boost::math::tools::promote_args<T1__>::type
36: integrate(const int& n,
37: const T1__& z, std::ostream* pstream__);
38:
39:
40: #include “C:/Users/nhatlth/Documents/exp_int.hpp”
41: class model29dc489d44f4_external : public prob_grad {
42: private:
43: public:
44: model29dc489d44f4_external(stan::io::var_context& context__,
45: std::ostream* pstream__ = 0)
46: : prob_grad(0) {
47: ctor_body(context__, 0, pstream__);
48: }
49:
50: model29dc489d44f4_external(stan::io::var_context& context__,
51: unsigned int random_seed__,
52: std::ostream* pstream__ = 0)
53: : prob_grad(0) {
54: ctor_body(context__, random_seed__, pstream__);
55: }
56:
57: void ctor_body(stan::io::var_context& context__,
58: unsigned int random_seed__,
59: std::ostream* pstream__) {
60: boost::ecuyer1988 base_rng__ =
61: stan::services::util::create_rng(random_seed__, 0);
62: (void) base_rng__; // suppress unused var warning
63:
64: current_statement_begin__ = -1;
65:
66: static const char* function__ = “model29dc489d44f4_external_namespace::model29dc489d44f4_external”;
67: (void) function__; // dummy to suppress unused var warning
68: size_t pos__;
69: (void) pos__; // dummy to suppress unused var warning
70: std::vector vals_i__;
71: std::vector vals_r__;
72: double DUMMY_VAR__(std::numeric_limits::quiet_NaN());
73: (void) DUMMY_VAR__; // suppress unused var warning
74:
75: // initialize member variables
76: try {
77:
78: // validate, data variables
79: // initialize data variables
80:
81:
82: // validate transformed data
83:
84: // validate, set parameter ranges
85: num_params_r__ = 0U;
86: param_ranges_i__.clear();
87: current_statement_begin__ = 5;
88: ++num_params_r__;
89: } catch (const std::exception& e) {
90: stan::lang::rethrow_located(e, current_statement_begin__, prog_reader__());
91: // Next line prevents compiler griping about no return
92: throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ");
93: }
94: }
95:
96: ~model29dc489d44f4_external() { }
97:
98:
99: void transform_inits(const stan::io::var_context& context__,
100: std::vector& params_i__,
101: std::vector& params_r__,
102: std::ostream pstream__) const {
103: stan::io::writer writer__(params_r__,params_i__);
104: size_t pos__;
105: (void) pos__; // dummy call to supress warning
106: std::vector vals_r__;
107: std::vector vals_i__;
108:
109: if (!(context__.contains_r(“x”)))
110: throw std::runtime_error(“variable x missing”);
111: vals_r__ = context__.vals_r(“x”);
112: pos__ = 0U;
113: context__.validate_dims(“initialization”, “x”, “double”, context__.to_vec());
114: double x(0);
115: x = vals_r__[pos__++];
116: try {
117: writer__.scalar_unconstrain(x);
118: } catch (const std::exception& e) {
119: throw std::runtime_error(std::string("Error transforming variable x: ") + e.what());
120: }
121:
122: params_r__ = writer__.data_r();
123: params_i__ = writer__.data_i();
124: }
125:
126: void transform_inits(const stan::io::var_context& context,
127: Eigen::Matrix<double,Eigen::Dynamic,1>& params_r,
128: std::ostream pstream__) const {
129: std::vector params_r_vec;
130: std::vector params_i_vec;
131: transform_inits(context, params_i_vec, params_r_vec, pstream__);
132: params_r.resize(params_r_vec.size());
133: for (int i = 0; i < params_r.size(); ++i)
134: params_r(i) = params_r_vec[i];
135: }
136:
137:
138: template <bool propto__, bool jacobian__, typename T__>
139: T__ log_prob(vector<T__>& params_r__,
140: vector& params_i__,
141: std::ostream pstream__ = 0) const {
142:
143: T__ DUMMY_VAR__(std::numeric_limits::quiet_NaN());
144: (void) DUMMY_VAR__; // suppress unused var warning
145:
146: T__ lp__(0.0);
147: stan::math::accumulator<T__> lp_accum__;
148:
149: try {
150: // model parameters
151: stan::io::reader<T__> in__(params_r__,params_i__);
152:
153: T__ x;
154: (void) x; // dummy to suppress unused var warning
155: if (jacobian__)
156: x = in__.scalar_constrain(lp__);
157: else
158: x = in__.scalar_constrain();
159:
160:
161: // transformed parameters
162:
163:
164:
165: // validate transformed parameters
166:
167: const char* function__ = “validate transformed params”;
168: (void) function__; // dummy to suppress unused var warning
169:
170: // model body
171:
172: current_statement_begin__ = 6;
173: lp_accum__.add(normal_log<propto__>(x, 0, 1));
174:
175: } catch (const std::exception& e) {
176: stan::lang::rethrow_located(e, current_statement_begin__, prog_reader__());
177: // Next line prevents compiler griping about no return
178: throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ");
179: }
180:
181: lp_accum__.add(lp__);
182: return lp_accum__.sum();
183:
184: } // log_prob()
185:
186: template <bool propto, bool jacobian, typename T_>
187: T_ log_prob(Eigen::Matrix<T_,Eigen::Dynamic,1>& params_r,
188: std::ostream pstream = 0) const {
189: std::vector<T_> vec_params_r;
190: vec_params_r.reserve(params_r.size());
191: for (int i = 0; i < params_r.size(); ++i)
192: vec_params_r.push_back(params_r(i));
193: std::vector vec_params_i;
194: return log_prob<propto,jacobian,T_>(vec_params_r, vec_params_i, pstream);
195: }
196:
197:
198: void get_param_names(std::vectorstd::string& names__) const {
199: names__.resize(0);
200: names__.push_back(“x”);
201: names__.push_back(“y”);
202: }
203:
204:
205: void get_dims(std::vector<std::vector<size_t> >& dimss__) const {
206: dimss__.resize(0);
207: std::vector<size_t> dims__;
208: dims__.resize(0);
209: dimss__.push_back(dims__);
210: dims__.resize(0);
211: dimss__.push_back(dims__);
212: }
213:
214: template
215: void write_array(RNG& base_rng__,
216: std::vector& params_r__,
217: std::vector& params_i__,
218: std::vector& vars__,
219: bool include_tparams__ = true,
220: bool include_gqs__ = true,
221: std::ostream pstream__ = 0) const {
222: vars__.resize(0);
223: stan::io::reader in__(params_r__,params_i__);
224: static const char function__ = “model29dc489d44f4_external_namespace::write_array”;
225: (void) function__; // dummy to suppress unused var warning
226: // read-transform, write parameters
227: double x = in__.scalar_constrain();
228: vars__.push_back(x);
229:
230: if (!include_tparams__) return;
231: // declare and define transformed parameters
232: double lp__ = 0.0;
233: (void) lp__; // dummy to suppress unused var warning
234: stan::math::accumulator lp_accum__;
235:
236: double DUMMY_VAR__(std::numeric_limits::quiet_NaN());
237: (void) DUMMY_VAR__; // suppress unused var warning
238:
239: try {
240:
241:
242:
243: // validate transformed parameters
244:
245: // write transformed parameters
246:
247: if (!include_gqs__) return;
248: // declare and define generated quantities
249: current_statement_begin__ = 8;
250: double y(0.0);
251: (void) y; // dummy to suppress unused var warning
252:
253: stan::math::initialize(y, std::numeric_limits::quiet_NaN());
254: stan::math::fill(y,DUMMY_VAR__);
255:
256:
257: current_statement_begin__ = 9;
258: stan::math::assign(y, (x + integrate(0,-(1.0), pstream__)));
259:
260: // validate generated quantities
261: current_statement_begin__ = 8;
262:
263: // write generated quantities
264: vars__.push_back(y);
265:
266: } catch (const std::exception& e) {
267: stan::lang::rethrow_located(e, current_statement_begin__, prog_reader__());
268: // Next line prevents compiler griping about no return
269: throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ");
270: }
271: }
272:
273: template
274: void write_array(RNG& base_rng,
275: Eigen::Matrix<double,Eigen::Dynamic,1>& params_r,
276: Eigen::Matrix<double,Eigen::Dynamic,1>& vars,
277: bool include_tparams = true,
278: bool include_gqs = true,
279: std::ostream pstream = 0) const {
280: std::vector params_r_vec(params_r.size());
281: for (int i = 0; i < params_r.size(); ++i)
282: params_r_vec[i] = params_r(i);
283: std::vector vars_vec;
284: std::vector params_i_vec;
285: write_array(base_rng,params_r_vec,params_i_vec,vars_vec,include_tparams,include_gqs,pstream);
286: vars.resize(vars_vec.size());
287: for (int i = 0; i < vars.size(); ++i)
288: vars(i) = vars_vec[i];
289: }
290:
291: static std::string model_name() {
292: return “model29dc489d44f4_external”;
293: }
294:
295:
296: void constrained_param_names(std::vectorstd::string& param_names__,
297: bool include_tparams__ = true,
298: bool include_gqs__ = true) const {
299: std::stringstream param_name_stream__;
300: param_name_stream__.str(std::string());
301: param_name_stream__ << “x”;
302: param_names__.push_back(param_name_stream__.str());
303:
304: if (!include_gqs__ && !include_tparams__) return;
305:
306: if (!include_gqs__) return;
307: param_name_stream__.str(std::string());
308: param_name_stream__ << “y”;
309: param_names__.push_back(param_name_stream__.str());
310: }
311:
312:
313: void unconstrained_param_names(std::vectorstd::string& param_names__,
314: bool include_tparams__ = true,
315: bool include_gqs__ = true) const {
316: std::stringstream param_name_stream__;
317: param_name_stream__.str(std::string());
318: param_name_stream__ << “x”;
319: param_names__.push_back(param_name_stream__.str());
320:
321: if (!include_gqs__ && !include_tparams__) return;
322:
323: if (!include_gqs__) return;
324: param_name_stream__.str(std::string());
325: param_name_stream__ << “y”;
326: param_names__.push_back(param_name_stream__.str());
327: }
328:
329: }; // model
330:
331: }
332:
333: typedef model29dc489d44f4_external_namespace::model29dc489d44f4_external stan_model;
334:
335: #include <rstan/rstaninc.hpp>
336: /
337: * Define Rcpp Module to expose stan_fit’s functions to R.
338: */
339: RCPP_MODULE(stan_fit4model29dc489d44f4_external_mod){
340: Rcpp::class_<rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external,
341: boost::random::ecuyer1988> >(“stan_fit4model29dc489d44f4_external”)
342: // .constructorRcpp::List()
343: .constructor<SEXP, SEXP, SEXP>()
344: // .constructor<SEXP, SEXP>()
345: .method(“call_sampler”,
346: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::call_sampler)
347: .method(“param_names”,
348: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::param_names)
349: .method(“param_names_oi”,
350: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::param_names_oi)
351: .method(“param_fnames_oi”,
352: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::param_fnames_oi)
353: .method(“param_dims”,
354: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::param_dims)
355: .method(“param_dims_oi”,
356: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::param_dims_oi)
357: .method(“update_param_oi”,
358: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::update_param_oi)
359: .method(“param_oi_tidx”,
360: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::param_oi_tidx)
361: .method(“grad_log_prob”,
362: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::grad_log_prob)
363: .method(“log_prob”,
364: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::log_prob)
365: .method(“unconstrain_pars”,
366: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::unconstrain_pars)
367: .method(“constrain_pars”,
368: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::constrain_pars)
369: .method(“num_pars_unconstrained”,
370: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::num_pars_unconstrained)
371: .method(“unconstrained_param_names”,
372: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::unconstrained_param_names)
373: .method(“constrained_param_names”,
374: &rstan::stan_fit<model29dc489d44f4_external_namespace::model29dc489d44f4_external, boost::random::ecuyer1988>::constrained_param_names)
375: ;
376: }
377:
378: // declarations
379: extern “C” {
380: SEXP file29dc94a1a87( ) ;
381: }
382:
383: // definition
384:
385: SEXP file29dc94a1a87( ){
386: return Rcpp::wrap(“external”);
387: }
388:
389:
Error in compileCode(f, code, language = language, verbose = verbose) :
Compilation ERROR, function(s)/method(s) not created! In file included from C:/Users/nhatlth/Documents/R/win-library/3.3/BH/include/boost/config.hpp:39:0,
from C:/Users/nhatlth/Documents/R/win-library/3.3/BH/include/boost/math/tools/config.hpp:13,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/core/var.hpp:7,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/core/gevv_vvv_vari.hpp:5,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/core.hpp:12,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math/rev/mat.hpp:4,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/stan/math.hpp:4,
from C:/Users/nhatlth/Documents/R/win-library/3.3/StanHeaders/include/src/stan/model/model_header.hpp:4,
from file29d
I don’t know what wrong with this? It has the similar formular as cyl_bessel_k function. Does anyone know how to solve this issue?

Thank you very much!
Nhat

25

Can you quote the instructions from that vignette that say that?

26

Hi all,

I have a solution. I put it here.
mc2 <-
’
functions {
//Adaptive Quadrature for simpson rule

real simp_adapt(real a,real b,real nu);
real adsimp(real a,real b,real nu,real fa,real fm,real fb,real V0,real delta);
real f(real x,real nu);
real f(real x,real nu){
real Q=exp(xnu);
return Q;
}
real simp_adapt(real a,real b,real nu){
if((anu>700)||(bnu>700)){
return exp(700);
}else{
real fa = f(a,nu);
real fm = f((a+b)/2,nu);
real fb = f(b,nu);
real V0 = (fa + 4fm + fb)(b-a)/6;
real delta;// error tolerance fixed
delta=1e-6;
return adsimp(a,b,nu,fa,fm,fb,V0,delta);
}
}
real adsimp(real a,real b,real nu,real fa,real fm,real fb,real V0,real delta){
real h = b-a;
real f1 = f(a + h/4,nu); // left-half midpoint
real f2 = f(b - h/4,nu); // right-half midpoint
real sl = h(fa + 4f1 + fm)/12 ;// left-half estimate
real sr = h(fm + 4*f2 + fb)/12;// right-half estimate
real s2 = sl+sr;
real err = (s2-V0)/15;
real V;
real m;
real V1;
real V2;
if (fabs(err)<delta) {
// estimate is within tolerance, so accept it
V = s2 + err;
}else{ //% split interval into two pieces
m = a + h/2;
V1 = adsimp(a,m,nu,fa,f1,fm,sl,delta/2);
V2 = adsimp(m,b,nu,fm,f2,fb,sr,delta/2);
V = V1 + V2;
};
return V;
}
}
parameters{
real x;
}

model {

} // use the sim_adapt() function somehow

’

model_tmp2<-stan_model(model_code = mc2)

Cheers,
Nhat