Added compile function.
authorJens Vollinga <vollinga@thep.physik.uni-mainz.de>
Wed, 7 Dec 2005 18:19:44 +0000 (18:19 +0000)
committerJens Vollinga <vollinga@thep.physik.uni-mainz.de>
Wed, 7 Dec 2005 18:19:44 +0000 (18:19 +0000)
configure.ac
doc/examples/compile1.cpp [new file with mode: 0644]
doc/examples/ginac-examples.texi
ginac/Makefile.am
ginac/excompiler.cpp [new file with mode: 0644]
ginac/excompiler.h [new file with mode: 0644]
ginac/ginac.h
tools/Makefile.am
tools/ginac-excompiler.in [new file with mode: 0644]

index 7005c94..47dcfdd 100644 (file)
@@ -116,6 +116,16 @@ if test "x$CONFIG_RUSAGE" = "xno"; then
     AC_CHECK_HEADER(ctime, , GINAC_ERROR([The standard <ctime> header file could not be found.]))
 fi
 
+dnl Check for dl library (needed for GiNaC::compile).
+AC_CHECK_LIB(dl, dlopen,
+       [
+               DL_LIBS="-ldl"
+               AC_SUBST(DL_LIBS)
+               AC_DEFINE(HAVE_LIBDL, 1, [set to 1 if you have a working libdl installed.])
+       ], 
+       GINAC_WARNING([libdl not found. GiNaC::compile will be disabled.]))
+LIBS="$LIBS $DL_LIBS"
+
 dnl We need to have Bruno Haible's CLN installed.
 dnl (CLN versions >= 1.1.0 must have installed cln.m4 at a visible place,
 dnl which provides this macro):
@@ -154,6 +164,7 @@ ginsh/Makefile
 ginsh/ginsh.1
 tools/Makefile
 tools/viewgar.1
+tools/ginac-excompiler
 doc/Makefile
 doc/examples/Makefile
 doc/tutorial/Makefile
diff --git a/doc/examples/compile1.cpp b/doc/examples/compile1.cpp
new file mode 100644 (file)
index 0000000..e2fa8e0
--- /dev/null
@@ -0,0 +1,45 @@
+#include <iostream>
+using namespace std;
+#include <ginac/ginac.h>
+using namespace GiNaC;
+// Yes, we are using CUBA (should be installed on the system!)
+#include "cuba.h"
+
+int main()
+{
+       // Let the user enter a expression
+       symbol x("x"), y("y");
+       string s;
+       cout << "Enter an expression containing 'x' and/or 'y': ";
+       cin >> s;
+       // Expression now in expr
+       ex expr(s, lst(x,y));
+       cout << "start integration of " << expr << " ..." << endl;
+
+       // Some definitions for VEGAS 
+       #define NDIM 2
+       #define NCOMP 1
+       #define EPSREL 1e-3
+       #define EPSABS 1e-12
+       #define VERBOSE 0
+       #define MINEVAL 0
+       #define MAXEVAL 50000
+       #define NSTART 1000
+       #define NINCREASE 500
+
+       // Some variables for VEGAS
+       int comp, nregions, neval, fail;
+       double integral[NCOMP], error[NCOMP], prob[NCOMP];
+
+       // Starting VEGAS
+       // By invocation of compile() the expression in expr is converted into the
+       // appropriate function pointer
+       Vegas(NDIM, NCOMP, compile(lst(expr), lst(x,y)), EPSREL, EPSABS, VERBOSE,
+             MINEVAL, MAXEVAL, NSTART, NINCREASE, &neval, &fail, integral, error, prob);
+
+       // Show the result
+       cout << "result: " << integral[0] << endl;
+
+       return 0;
+}
index 200c81e..bd14dab 100644 (file)
@@ -30,4 +30,11 @@ This is a collection of code examples using GiNaC.
 
 Two expression are stored in an archive on the disk and are restored again.
 
+@chapter Monte Carlo Integration
+
+@section Using VEGAS from CUBA @uref{compile1.cpp, (source)}
+
+An expression in two variables is integrated by using VEGAS from the
+@uref{http://www.feynarts.de/cuba/, CUBA library}.
+
 @bye
index 0236f9b..7ae8455 100644 (file)
@@ -2,8 +2,8 @@
 
 lib_LTLIBRARIES = libginac.la
 libginac_la_SOURCES = add.cpp archive.cpp basic.cpp clifford.cpp color.cpp \
-  constant.cpp ex.cpp expair.cpp expairseq.cpp exprseq.cpp fail.cpp \
-  fderivative.cpp function.cpp idx.cpp indexed.cpp inifcns.cpp \
+  constant.cpp ex.cpp excompiler.cpp expair.cpp expairseq.cpp exprseq.cpp \
+  fail.cpp fderivative.cpp function.cpp idx.cpp indexed.cpp inifcns.cpp \
   inifcns_trans.cpp inifcns_gamma.cpp inifcns_nstdsums.cpp \
   integral.cpp lst.cpp matrix.cpp mul.cpp ncmul.cpp normal.cpp numeric.cpp \
   operators.cpp power.cpp registrar.cpp relational.cpp remember.cpp \
@@ -13,9 +13,9 @@ libginac_la_SOURCES = add.cpp archive.cpp basic.cpp clifford.cpp color.cpp \
 libginac_la_LDFLAGS = -version-info $(LT_VERSION_INFO) -release $(LT_RELEASE)
 ginacincludedir = $(includedir)/ginac
 ginacinclude_HEADERS = ginac.h add.h archive.h assertion.h basic.h class_info.h \
-  clifford.h color.h constant.h container.h ex.h expair.h expairseq.h exprseq.h \
-  fail.h fderivative.h flags.h function.h hash_map.h idx.h indexed.h inifcns.h \
-  integral.h lst.h matrix.h mul.h ncmul.h normal.h numeric.h operators.h \
+  clifford.h color.h constant.h container.h ex.h excompiler.h expair.h expairseq.h \
+  exprseq.h fail.h fderivative.h flags.h function.h hash_map.h idx.h indexed.h \
+  inifcns.h integral.h lst.h matrix.h mul.h ncmul.h normal.h numeric.h operators.h \
   power.h print.h pseries.h ptr.h registrar.h relational.h structure.h \
   symbol.h symmetry.h tensor.h tinfos.h version.h wildcard.h
 AM_LFLAGS = -Pginac_yy -olex.yy.c
diff --git a/ginac/excompiler.cpp b/ginac/excompiler.cpp
new file mode 100644 (file)
index 0000000..4f121d7
--- /dev/null
@@ -0,0 +1,170 @@
+/** @file excompiler.cpp
+ *
+ *  Class to facilitate the conversion of a ex to a function pointer suited for
+ *  fast numerical integration. */
+
+/*
+ *  GiNaC Copyright (C) 1999-2005 Johannes Gutenberg University Mainz, Germany
+ *
+ *  This program is free software; you can redistribute it and/or modify
+ *  it under the terms of the GNU General Public License as published by
+ *  the Free Software Foundation; either version 2 of the License, or
+ *  (at your option) any later version.
+ *
+ *  This program is distributed in the hope that it will be useful,
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ *  GNU General Public License for more details.
+ *
+ *  You should have received a copy of the GNU General Public License
+ *  along with this program; if not, write to the Free Software
+ *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
+ */
+
+#include "excompiler.h"
+
+#include <dlfcn.h>
+#include <stdexcept>
+#include <fstream>
+#include <sstream>
+#include <string>
+#include <vector>
+
+#include "ex.h"
+#include "lst.h"
+#include "operators.h"
+#include "relational.h"
+#include "symbol.h"
+
+namespace GiNaC {
+
+#ifdef HAVE_LIBDL
+       
+class excompiler
+{
+       struct filedesc
+       {
+               void* module;
+               std::string name;
+       };
+public:
+       excompiler() {};
+       ~excompiler()
+       {
+               for (std::vector<filedesc>::const_iterator it = filelist.begin(); it != filelist.end(); ++it) {
+                       dlclose(it->module);
+                       std::string strremove = "rm " + it->name;
+                       system(strremove.c_str());
+               }
+       }
+       void add(void* module, std::string name)
+       {
+               filedesc fd;
+               fd.module = module;
+               fd.name = name;
+               filelist.push_back(fd);
+       }
+       std::vector<filedesc> filelist;
+};
+
+excompiler _exc;
+
+FP_dim1 compile(const ex& expr, const symbol& sym)
+{
+       symbol argx("argx");
+       ex expr_with_x = expr.subs(lst(sym==argx));
+
+       char filename[] = "/tmp/GiNaCXXXXXX";
+
+       int fno = mkstemp(filename);
+
+       std::ofstream ofs(filename);
+
+       ofs << "#include <stddef.h> " << std::endl;
+       ofs << "#include <stdlib.h> " << std::endl;
+       ofs << "#include <math.h> " << std::endl;
+       ofs << std::endl;
+       ofs << "double compiled_ex(double argx)" << std::endl;
+       ofs << "{" << std::endl;
+       ofs << "double res = ";
+       expr_with_x.print(GiNaC::print_csrc_double(ofs));
+       ofs << ";" << std::endl;
+       ofs << "return(res); " << std::endl;
+       ofs << "}" << std::endl;
+
+       ofs.close();
+
+       std::string strcompile = "ginac-excompiler " + std::string(filename);
+       system(strcompile.c_str());
+
+       std::string strremove = "rm " + std::string(filename) + " " + std::string(filename) + ".o";
+       system(strremove.c_str());
+
+       std::string strsofile = std::string(filename) + ".so";
+       void* module = NULL;
+       module = dlopen(strsofile.c_str(), RTLD_NOW);
+       if (module == NULL)     {
+               throw std::runtime_error("excompiler: could not open compiled module!");
+       }
+
+       _exc.add(module, strsofile);
+
+       return (FP_dim1) dlsym(module, "compiled_ex");
+}
+
+FP_cuba compile(const lst& exprs, const lst& syms)
+{
+       lst replacements;
+       for (int count=0; count<syms.nops(); ++count) {
+               std::ostringstream s;
+               s << "a[" << count << "]";
+               replacements.append(syms.op(count) == symbol(s.str()));
+       }
+
+       std::vector<ex> expr_with_cname;
+       for (int count=0; count<exprs.nops(); ++count) {
+               expr_with_cname.push_back(exprs.op(count).subs(replacements));
+       }
+
+       char filename[] = "/tmp/GiNaCXXXXXX";
+
+       int fno = mkstemp(filename);
+
+       std::ofstream ofs(filename);
+
+       ofs << "#include <stddef.h> " << std::endl;
+       ofs << "#include <stdlib.h> " << std::endl;
+       ofs << "#include <math.h> " << std::endl;
+       ofs << std::endl;
+       ofs << "void compiled_ex(const int* an, const double a[], const int* fn, double f[])" << std::endl;
+       ofs << "{" << std::endl;
+       for (int count=0; count<exprs.nops(); ++count) {
+               ofs << "f[" << count << "] = ";
+               expr_with_cname[count].print(GiNaC::print_csrc_double(ofs));
+               ofs << ";" << std::endl;
+       }
+       ofs << "}" << std::endl;
+
+       ofs.close();
+
+       std::string strcompile = "ginac-excompiler " + std::string(filename);
+       system(strcompile.c_str());
+
+       std::string strremove = "rm " + std::string(filename) + " " + std::string(filename) + ".o";
+       system(strremove.c_str());
+
+       std::string strsofile = std::string(filename) + ".so";
+       void* module = NULL;
+       module = dlopen(strsofile.c_str(), RTLD_NOW);
+       if (module == NULL)     {
+               throw std::runtime_error("excompiler: could not open compiled module!");
+       }
+
+       _exc.add(module, strsofile);
+
+       return (FP_cuba) dlsym(module, "compiled_ex");
+}
+
+#endif
+
+} // namespace GiNaC
diff --git a/ginac/excompiler.h b/ginac/excompiler.h
new file mode 100644 (file)
index 0000000..0c1ca3a
--- /dev/null
@@ -0,0 +1,47 @@
+/** @file excompiler.h
+ *
+ *  Class to facilitate the conversion of a ex to a function pointer suited for
+ *  fast numerical integration. */
+
+/*
+ *  GiNaC Copyright (C) 1999-2005 Johannes Gutenberg University Mainz, Germany
+ *
+ *  This program is free software; you can redistribute it and/or modify
+ *  it under the terms of the GNU General Public License as published by
+ *  the Free Software Foundation; either version 2 of the License, or
+ *  (at your option) any later version.
+ *
+ *  This program is distributed in the hope that it will be useful,
+ *  but WITHOUT ANY WARRANTY; without even the implied warranty of
+ *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ *  GNU General Public License for more details.
+ *
+ *  You should have received a copy of the GNU General Public License
+ *  along with this program; if not, write to the Free Software
+ *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
+ */
+
+#ifndef __GINAC_EXCOMPILER_H__
+#define __GINAC_EXCOMPILER_H__
+
+#include "basic.h"
+#include "config.h"
+#include "ex.h"
+
+namespace GiNaC {
+
+#ifdef HAVE_LIBDL
+
+typedef double (*FP_dim1) (double);
+
+FP_dim1 compile(const ex& expr, const symbol& sym);
+
+typedef void (*FP_cuba) (const int*, const double[], const int*, double[]);
+
+FP_cuba compile(const lst& exprs, const lst& syms);
+
+#endif
+
+} // namespace GiNaC
+
+#endif // ndef __GINAC_EXCOMPILER_H__
index 4f4f23a..0bee28e 100644 (file)
@@ -66,6 +66,8 @@
 #include "color.h"
 #include "clifford.h"
 
+#include "excompiler.h"
+
 #ifdef __MAKECINT__
 #pragma link C++ nestedclass;
 #pragma link C++ nestedtypedef;
index a058b74..52280a3 100644 (file)
@@ -4,6 +4,10 @@ bin_PROGRAMS = viewgar
 viewgar_SOURCES = viewgar.cpp
 viewgar_LDADD = ../ginac/libginac.la
 
+bin_SCRIPTS = ginac-excompiler
+
 INCLUDES = -I$(srcdir)/../ginac -I../ginac
 
 man_MANS = viewgar.1
+
+EXTRA_DIST = ginac-excompiler.in
diff --git a/tools/ginac-excompiler.in b/tools/ginac-excompiler.in
new file mode 100644 (file)
index 0000000..e547edd
--- /dev/null
@@ -0,0 +1,3 @@
+#!/bin/sh
+@CC@ -x c -fPIC -o $1.o -c $1
+@CC@ -shared -o $1.so $1.o