Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
Minimization.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_r
3/// \notebook -nodraw
4/// Example based in
5/// http://root.cern.ch/root/html/tutorials/fit/NumericalMinimization.C.html
6/// http://stat.ethz.ch/R-manual/R-devel/library/stats/html/optim.html
7///
8/// \macro_code
9///
10/// \author Omar Zapata
11
12#include<TRInterface.h>
13
14//in the next function the *double pointer must be changed by a TVectorD,
15//because the pointer has no meaning in R enviroment.
16Double_t RosenBrock(const TVectorD xx )
17{
18 const Double_t x = xx[0];
19 const Double_t y = xx[1];
20 const Double_t tmp1 = y-x*x;
21 const Double_t tmp2 = 1-x;
22 return 100*tmp1*tmp1+tmp2*tmp2;
23}
24
25TVectorD RosenBrockGrad(const TVectorD xx )
26{
27 const Double_t x = xx[0];
28 const Double_t y = xx[1];
29 TVectorD grad(2);
30 grad[0]=-400 * x * (y - x * x) - 2 * (1 - x);
31 grad[1]=200 * (y - x * x);
32 return grad;
33}
34
35void Minimization()
36{
38
39 //passsing RosenBrock function to R
40 r["RosenBrock"]=ROOT::R::TRFunctionExport(RosenBrock);
41
42 //passsing RosenBrockGrad function to R
43 r["RosenBrockGrad"]=ROOT::R::TRFunctionExport(RosenBrockGrad);
44 //the option "method" could be "Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN","Brent"
45
46 //the option "control" lets you put some constraints like
47 //"maxit" The maximum number of iterations.
48 //"abstol" The absolute convergence tolerance.
49 r.Execute("result <- optim( c(0.01,0.01), RosenBrock,method='BFGS',control = list(maxit = 1000000) )");
50 //"reltol" Relative convergence tolerance.
51
52 //Getting results from R
53 TVectorD min=r.Eval("result$par");
54
55 std::cout.precision(8);
56 //printing results
57 std::cout<<"-----------------------------------------"<<std::endl;
58 std::cout<<"Minimum x="<<min[0]<<" y="<<min[1]<<std::endl;
59 std::cout<<"Value at minimum ="<<RosenBrock(min)<<std::endl;
60
61 //using the gradient
62 r.Execute("optimHess(result$par, RosenBrock, RosenBrockGrad)");
63 r.Execute("hresult <- optim(c(-1.2,1), RosenBrock, NULL, method = 'BFGS', hessian = TRUE)");
64 //getting the min calculated with the gradient
65 TVectorD hmin=r.Eval("hresult$par");
66
67 //printing results
68 std::cout<<"-----------------------------------------"<<std::endl;
69 std::cout<<"Minimization with the Gradient"<<std::endl;
70 std::cout<<"Minimum x="<<hmin[0]<<" y="<<hmin[1]<<std::endl;
71 std::cout<<"Value at minimum ="<<RosenBrock(hmin)<<std::endl;
72}
ROOT::R::TRInterface & r
Definition Object.C:4
double Double_t
Definition RtypesCore.h:59
This is a class to pass functions from ROOT to R.
ROOT R was implemented using the R Project library and the modules Rcpp and RInside
void Execute(const TString &code)
Method to eval R code.
static TRInterface & Instance()
static method to get an TRInterface instance reference
Int_t Eval(const TString &code, TRObject &ans)
Method to eval R code and you get the result in a reference to TRObject.
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17