Brief Update, Avraham's advice helped me get a better idea of where to start. What I am now trying to do is learn from this post, which explains how to use external libraries in RCPP ( http://stackoverflow.com/questions/13995266/using-3rd-party-header-files-with-rcpp). And secondly, try to use that to implement this optimize function ( http://www.boost.org/doc/libs/1_57_0/libs/math/doc/html/math_toolkit/internals1/minima.html ).
If anyone has additional links or input I would be grateful to receive it, I do not expect anyone to hold my hand through this. Thank you again for your time and help, Simon On Mon, Dec 22, 2014 at 9:25 PM, Avraham Adler <avraham.ad...@gmail.com> wrote: > Hello, Simon. > > I ran into a similar problem before (with uniroot > <https://stat.ethz.ch/R-manual/R-devel/library/stats/html/uniroot.html>), > and what I found was that the function in R is actually written in C. I > found the original code, but not being a real programmer, I stopped working > on what I was doing and put it off until the point where I have enough time > to develop the skill to port that into RCpp (maybe 2047?). In your case, > optimize() > <https://stat.ethz.ch/R-manual/R-devel/library/stats/html/optimize.html>itself > is a translation of an algorithm by Brent in FORTRAN which can be found here > (textfile) <http://www.netlib.org/fmm/fmin.f>. Perhaps you can translate > it into C++ and call it directly? > > Avi > > On Tue, Dec 23, 2014 at 12:04 AM, Simon Riddell <simon...@gmail.com> > wrote: > >> Hello, >> >> I have been judiciously using RCPP for six months, and have had my >> numerous questions answered so far from previous threads. I have now run >> into an issue I cannot get a handle on: >> >> I have converted a fairly large (recursive) Kalman filter function from R >> to C++ using RCPP. This function is then called within the R optim() >> function, subject to constraints, which estimates a set of ~30 parameters. >> >> My issue occurs within the C++ Kalman Filter function, where, I use the >> optimize() function to calculate the yield to maturity rate of various >> bonds. I do this by calling back to the R environment to access the >> optimize() function. Below is the R Function I create to be used within the >> Kalman filter, and below this R function is my method for calling it within >> the C++ code. To complicate matters further, the R Function calls a C++ >> Function. To clarify: The Kalman Filter C++ code calls an R function, and >> this R Function calls an additional separate C++ function. (Code included >> below) >> >> As I iterate the Kalman filter it runs perfectly for anywhere from 30 >> minutes to six hours (and produces correct output when matched to the R >> code), but it inevitably crashes. From past reading I have done, Dirk has >> before mentioned that calling an R function many times within C++ is >> usually not a good idea. As a result I suspect this is my issue (the error >> codes vary, sometimes mentioning numerical errors, other times recursive >> errors, other times random RCPP error codes -- I can document and provide >> them if needed) >> >> My biggest impasse is I cannot figure out a way to complete this without >> calling the R optimize() function, and cannot find any RCPP optimize >> functions to use instead. Thank you for reading. >> >> >> *R Function (**IRNPV.CPP is the C++ function, which R optimizes the rate >> parameter over):* >> >> optimcpp<-function(CoupDate=1,coupNo=1,coup=1,price=1,rate=1) >> { >> m<-optimize(function(CoupDate,coupNo,coup,price,rate) *IRNPV.CPP* >> (CoupDate=CoupDate,coupNo=coupNo,coup=coup,*rate* >> )-price)^2,c(-0.05,0.2),tol=1e-20,CoupDate=CoupDate,coupNo=coupNo,coup=coup,price=price) >> m$minimum >> } >> >> *Accessing the R environment within the C++ code:* >> CPP.SRC <- ' >> using namespace arma; >> Rcpp::Environment global = Rcpp::Environment::global_env(); >> Function optimizecpp = global["optimcpp"]; >> //Various matrix computations to calculate CD, CN, Co, & Pricex[0] >> optimvec0 = optimizecpp(CD,CN,Co,pricex[0],Ra); >> ' >> >> *IRNPV.CPP Function *(What the R Optimize() function optimizes 'rate' >> over -- *Very* *Likely Unnecessary for Purposes of this Question)* >> >> IR.NPV.TIPS.CBF.SRC <- ' >> using namespace arma; >> >> double CD = Rcpp::as<double>(CoupDate); >> double CN = Rcpp::as<double>(coupNo); >> double RN = Rcpp::as<double>(rate); >> double Co = Rcpp::as<double>(coup); >> >> Rcpp::NumericVector LM; >> mat DiscountFunc; >> double length; >> >> double price; >> >> if (CN > 1) { >> >> length = floor((((CD+0.5*(CN-1))-CD)/0.50))+1; >> >> for (double i = CD; i <= (CD+0.5*(CN-1))+.05; i += 0.5) { >> LM.insert(LM.end(),i); >> } >> >> DiscountFunc.set_size(LM.size(), 1); >> DiscountFunc.fill(0); >> >> >> double k = 0; >> mat::row_iterator q = DiscountFunc.begin_row(0); >> mat::row_iterator w = DiscountFunc.end_row((LM.size()-1)); >> for(mat::row_iterator i=q; i!=w; ++i) >> { >> (*i) = exp(-RN*LM[k]); >> k = k + 1; >> } >> >> price = CD*Co*DiscountFunc[0]; >> >> for (int i=1; i<(LM.size()); ++i) { >> price = price+0.5*Co*DiscountFunc[i]; >> } >> >> >> } >> else { >> double DiscountFunc; >> DiscountFunc = exp(-RN*CD); >> price = (1+CD*Co)*DiscountFunc; >> } >> return Rcpp::wrap(price); >> >> ' >> >> >> IRNPV.CPP <- cxxfunction(signature(CoupDate="NumericVector", >> coupNo="NumericVector", coup="NumericVector",rate="NumericVector"), >> IR.NPV.TIPS.CBF.SRC, include=CMATH, plugin="RcppArmadillo") >> >> >> >> Thank you, >> Simon >> >> >> -- >> Simon Alexander Riddell >> Economic Research RA >> Federal Reserve Bank >> >> _______________________________________________ >> Rcpp-devel mailing list >> Rcpp-devel@lists.r-forge.r-project.org >> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel >> > > -- Simon Alexander Riddell London School of Economics linkedin.com/in/simonriddell <http://uk.linkedin.com/in/simonriddell>
_______________________________________________ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
