Bug fix in the routine H_evalf (a minus sign should not be forgotten).
[ginac.git] / doc / examples / compile2.cpp
1 #include <iostream>
2 #ifdef IN_GINAC
3 #include "ginac.h"
4 #else
5 #include <ginac/ginac.h>
6 #endif
7 // Yes, we are using CUBA (should be installed on the system!)
8 #include <cuba.h>
9
10 using namespace std;
11 using namespace GiNaC;
12
13 /*
14  * Demonstrates the use of compile_ex with the CUBA library.
15  *
16  * The user can enter an expression on the command line. This expression is
17  * compiled via compile_ex and integrated over the region 0 <= x,y <= 1 with
18  * the help of the CUBA library (http://www.feynarts.de/cuba).
19  *
20  */
21
22 int main()
23 {
24         // Let the user enter a expression
25         symbol x("x"), y("y");
26         string s;
27         cout << "Enter an expression containing 'x' and/or 'y': ";
28         cin >> s;
29         // Expression now in expr
30         ex expr(s, lst{x,y});
31  
32         cout << "start integration of " << expr << " ..." << endl;
33
34         // Some definitions for VEGAS 
35         #define NDIM 2
36         #define NCOMP 1
37         #define EPSREL 1e-3
38         #define EPSABS 1e-12
39         #define VERBOSE 0
40         #define MINEVAL 0
41         #define MAXEVAL 50000
42         #define NSTART 1000
43         #define NINCREASE 500
44
45         // Some variables for VEGAS
46         int comp, nregions, neval, fail;
47         double integral[NCOMP], error[NCOMP], prob[NCOMP];
48
49         // Our function pointer that points to the compiled ex
50         FUNCP_CUBA fp;
51
52         // Optionally, compile with custom compiler flags:
53         // setenv("CXXFLAGS", "-O3 -fomit-frame-pointer -ffast-math", 1);
54         compile_ex(lst{expr}, lst{x,y}, fp);
55
56         // Starting VEGAS
57         // By invocation of compile() the expression in expr is converted into the
58         // appropriate function pointer
59         Vegas(NDIM, NCOMP, fp,
60               EPSREL, EPSABS, VERBOSE, MINEVAL, MAXEVAL, NSTART, NINCREASE,
61               &neval, &fail, integral, error, prob);
62
63         // Show the result
64         cout << "result: " << integral[0] << endl;
65
66         return 0;
67 }