mslang.cpp 5.27 KB
Newer Older
1
#include "Algorithms/Vektor.h"
2
#include "Utilities/MSLang.h"
kraus's avatar
kraus committed
3

4 5
#include "Ippl.h"
#include "Utility/IpplTimings.h"
6

7 8
#include <gsl/gsl_rng.h>

kraus's avatar
kraus committed
9 10
#include <boost/regex.hpp>

11 12
#include <iostream>
#include <string>
13
#include <fstream>
14
#include <sstream>
15
#include <streambuf>
16
#include <cstdlib>
17
#include <cmath>
18
#include <set>
19

20
Ippl *ippl;
21
int main(int argc, char *argv[])
22
{
23 24 25 26
    ippl = new Ippl(argc, argv);

    if (argc < 4) {
        std::cout << "please provide the name of the file that contains your code, the method and the number of particles" << std::endl;
27 28
        return 1;
    }
29
    mslang::Function *fun;
30 31 32 33 34

    std::ifstream t(argv[1]);
    std::string str((std::istreambuf_iterator<char>(t)),
                    std::istreambuf_iterator<char>());

35 36 37
    const unsigned int method = (atoi(argv[2]) == 0 ? 0: 1);
    const unsigned int N = atoi(argv[3]);

38
    // std::string str("repeat( translate(union(rectangle(0.1, 0.1), ellipse(0.1, 0.1)), -0.01, -0.02), 2, 0.1, 0.2)");
39

kraus's avatar
kraus committed
40 41 42 43
    str = boost::regex_replace(str, boost::regex("//.*?\\n"), std::string(""), boost::match_default | boost::format_all);
    str = boost::regex_replace(str, boost::regex("\\s"), std::string(""), boost::match_default | boost::format_all);

     if (parse(str, fun)) {
44 45
        // fun->print(0);
        // std::cout << "\n" << std::endl;
kraus's avatar
kraus committed
46

47
        std::vector<mslang::Base*> baseBlocks;
kraus's avatar
kraus committed
48 49
        fun->apply(baseBlocks);

50
        std::ofstream out("test.gpl");
51 52
        for (mslang::Base* bfun: baseBlocks) {
            // bfun->print(0);
53
            bfun->computeBoundingBox();
54
            bfun->writeGnuplot(out);
kraus's avatar
kraus committed
55
        }
56
        out.close();
57 58 59

        if (baseBlocks.size() > 0) {
            Vector_t llc, urc;
60 61
            mslang::Base* first = baseBlocks.front();
            const mslang::BoundingBox &bb = first->bb_m;
62 63 64 65 66 67 68
            llc = Vector_t(bb.center_m[0] - 0.5 * bb.width_m,
                           bb.center_m[1] - 0.5 * bb.height_m,
                           0.0);
            urc = Vector_t(bb.center_m[0] + 0.5 * bb.width_m,
                           bb.center_m[1] + 0.5 * bb.height_m,
                           0.0);

69

70
            for (unsigned int i = 1; i < baseBlocks.size(); ++ i) {
71
                const mslang::BoundingBox &bb = baseBlocks[i]->bb_m;
72 73 74 75 76 77
                llc[0] = std::min(llc[0], bb.center_m[0] - 0.5 * bb.width_m);
                llc[1] = std::min(llc[1], bb.center_m[1] - 0.5 * bb.height_m);
                urc[0] = std::max(urc[0], bb.center_m[0] + 0.5 * bb.width_m);
                urc[1] = std::max(urc[1], bb.center_m[1] + 0.5 * bb.height_m);
            }

kraus's avatar
kraus committed
78 79 80 81 82 83 84
            double width = urc[0] - llc[0];
            double height = urc[1] - llc[1];
            llc[0] -= 1e-3 * width;
            urc[0] += 1e-3 * width;
            llc[1] -= 1e-3 * height;
            urc[1] += 1e-3 * width;

85 86
            mslang::QuadTree tree;
            tree.bb_m = mslang::BoundingBox(llc, urc);
87

88 89
            tree.objects_m.insert(tree.objects_m.end(), baseBlocks.begin(), baseBlocks.end());
            tree.buildUp();
90

91 92
            out.open("quadtree.gpl");
            tree.writeGnuplot(out);
93 94 95

            out.close();

96
            if (true) {
97 98
                out.open("particles.gpl");
                gsl_rng *rng = gsl_rng_alloc(gsl_rng_default);
99

100 101 102
                unsigned int n = 0;
                IpplTimings::TimerRef mainTimer = IpplTimings::getTimer("mainTimer");
                IpplTimings::startTimer(mainTimer);
103

104 105 106 107
                for (unsigned int i = 0; i < N; ++ i) {
                    Vector_t X(0.0);
                    X[0] = llc[0] + (urc[0] - llc[0]) * gsl_rng_uniform(rng);
                    X[1] = llc[1] + (urc[1] - llc[1]) * gsl_rng_uniform(rng);
108

109

110 111
                    if (method == 0) {
                        if (tree.isInside(X)) {
112 113 114 115
                            ++ n;
                            out << std::setw(14) << X[0]
                                << std::setw(14) << X[1]
                                << std::endl;
116 117 118 119 120 121 122 123 124 125
                        }
                    } else {
                        for (mslang::Base* func: baseBlocks) {
                            if (func->isInside(X)) {
                                ++ n;
                                out << std::setw(14) << X[0]
                                    << std::setw(14) << X[1]
                                    << std::endl;
                                break;
                            }
126
                        }
127 128
                    }
                }
129 130 131 132 133 134 135 136 137 138 139 140
                IpplTimings::stopTimer(mainTimer);

                std::cout << (double)n / N * 100 << " % of particles passed" << std::endl;

                std::map<std::string, unsigned int> characteristicValues;
                characteristicValues.insert(std::make_pair("method", method));
                characteristicValues.insert(std::make_pair("num particles", N));
                characteristicValues.insert(std::make_pair("num base functions", baseBlocks.size()));
                std::stringstream ss;
                ss << "timing__m_=_" << method << "__np_=_" << N << "__nbf_=_" << baseBlocks.size() << ".dat";
                IpplTimings::print(ss.str(), characteristicValues);
                gsl_rng_free(rng);
141
            }
kraus's avatar
kraus committed
142

143
            for (mslang::Base* func: baseBlocks) {
144 145
                delete func;
            }
kraus's avatar
kraus committed
146
        }
147 148

        out.close();
149
    }
kraus's avatar
kraus committed
150 151 152

    delete fun;

153 154
    return 0;
}