Improved algebraic substitutions.
[ginac.git] / ginac / excompiler.cpp
1 /** @file excompiler.cpp
2  *
3  *  Class to facilitate the conversion of a ex to a function pointer suited for
4  *  fast numerical integration. */
5
6 /*
7  *  GiNaC Copyright (C) 1999-2005 Johannes Gutenberg University Mainz, Germany
8  *
9  *  This program is free software; you can redistribute it and/or modify
10  *  it under the terms of the GNU General Public License as published by
11  *  the Free Software Foundation; either version 2 of the License, or
12  *  (at your option) any later version.
13  *
14  *  This program is distributed in the hope that it will be useful,
15  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17  *  GNU General Public License for more details.
18  *
19  *  You should have received a copy of the GNU General Public License
20  *  along with this program; if not, write to the Free Software
21  *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
22  */
23
24 #include "excompiler.h"
25
26 #include <dlfcn.h>
27 #include <stdexcept>
28 #include <fstream>
29 #include <sstream>
30 #include <string>
31 #include <vector>
32
33 #include "ex.h"
34 #include "lst.h"
35 #include "operators.h"
36 #include "relational.h"
37 #include "symbol.h"
38
39 namespace GiNaC {
40
41 #ifdef HAVE_LIBDL
42         
43 class excompiler
44 {
45         struct filedesc
46         {
47                 void* module;
48                 std::string name;
49         };
50 public:
51         excompiler() {};
52         ~excompiler()
53         {
54                 for (std::vector<filedesc>::const_iterator it = filelist.begin(); it != filelist.end(); ++it) {
55                         dlclose(it->module);
56                         std::string strremove = "rm " + it->name;
57                         system(strremove.c_str());
58                 }
59         }
60         void add(void* module, std::string name)
61         {
62                 filedesc fd;
63                 fd.module = module;
64                 fd.name = name;
65                 filelist.push_back(fd);
66         }
67         std::vector<filedesc> filelist;
68 };
69
70 excompiler _exc;
71
72 FP_dim1 compile(const ex& expr, const symbol& sym)
73 {
74         symbol argx("argx");
75         ex expr_with_x = expr.subs(lst(sym==argx));
76
77         char filename[] = "/tmp/GiNaCXXXXXX";
78
79         int fno = mkstemp(filename);
80
81         std::ofstream ofs(filename);
82
83         ofs << "#include <stddef.h> " << std::endl;
84         ofs << "#include <stdlib.h> " << std::endl;
85         ofs << "#include <math.h> " << std::endl;
86         ofs << std::endl;
87         ofs << "double compiled_ex(double argx)" << std::endl;
88         ofs << "{" << std::endl;
89         ofs << "double res = ";
90         expr_with_x.print(GiNaC::print_csrc_double(ofs));
91         ofs << ";" << std::endl;
92         ofs << "return(res); " << std::endl;
93         ofs << "}" << std::endl;
94
95         ofs.close();
96
97         std::string strcompile = "ginac-excompiler " + std::string(filename);
98         system(strcompile.c_str());
99
100         std::string strremove = "rm " + std::string(filename) + " " + std::string(filename) + ".o";
101         system(strremove.c_str());
102
103         std::string strsofile = std::string(filename) + ".so";
104         void* module = NULL;
105         module = dlopen(strsofile.c_str(), RTLD_NOW);
106         if (module == NULL)     {
107                 throw std::runtime_error("excompiler: could not open compiled module!");
108         }
109
110         _exc.add(module, strsofile);
111
112         return (FP_dim1) dlsym(module, "compiled_ex");
113 }
114
115 FP_cuba compile(const lst& exprs, const lst& syms)
116 {
117         lst replacements;
118         for (int count=0; count<syms.nops(); ++count) {
119                 std::ostringstream s;
120                 s << "a[" << count << "]";
121                 replacements.append(syms.op(count) == symbol(s.str()));
122         }
123
124         std::vector<ex> expr_with_cname;
125         for (int count=0; count<exprs.nops(); ++count) {
126                 expr_with_cname.push_back(exprs.op(count).subs(replacements));
127         }
128
129         char filename[] = "/tmp/GiNaCXXXXXX";
130
131         int fno = mkstemp(filename);
132
133         std::ofstream ofs(filename);
134
135         ofs << "#include <stddef.h> " << std::endl;
136         ofs << "#include <stdlib.h> " << std::endl;
137         ofs << "#include <math.h> " << std::endl;
138         ofs << std::endl;
139         ofs << "void compiled_ex(const int* an, const double a[], const int* fn, double f[])" << std::endl;
140         ofs << "{" << std::endl;
141         for (int count=0; count<exprs.nops(); ++count) {
142                 ofs << "f[" << count << "] = ";
143                 expr_with_cname[count].print(GiNaC::print_csrc_double(ofs));
144                 ofs << ";" << std::endl;
145         }
146         ofs << "}" << std::endl;
147
148         ofs.close();
149
150         std::string strcompile = "ginac-excompiler " + std::string(filename);
151         system(strcompile.c_str());
152
153         std::string strremove = "rm " + std::string(filename) + " " + std::string(filename) + ".o";
154         system(strremove.c_str());
155
156         std::string strsofile = std::string(filename) + ".so";
157         void* module = NULL;
158         module = dlopen(strsofile.c_str(), RTLD_NOW);
159         if (module == NULL)     {
160                 throw std::runtime_error("excompiler: could not open compiled module!");
161         }
162
163         _exc.add(module, strsofile);
164
165         return (FP_cuba) dlsym(module, "compiled_ex");
166 }
167
168 #endif
169
170 } // namespace GiNaC