Cpp_global_opt/Elbme.h at master · strazdas/Cpp_global_opt · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
#ifndef ELLEME_H
#define ELLEME_H
/* This algorithm was supposed to solve all inner problem, which consists of
 * - initializing problem dividing till criteria is met */
/* Efficient Lower bound Minimum Estimation ~ Elbme */
#include <algorithm>
#include <fstream>
#include <sstream>
#include <iostream>
#include <limits>
#include <math.h>
#include <string>
#include <sys/time.h>
#include <vector>
#include "utils.h"
#include "subsimplex.h"
// Depencencies are:   main  <  Disimplv  <  Simplex  <  Elbme  <  Subsimplex

using namespace std;

class Elbme {   // For efficiency do not store links to parents
    Elbme(const Elbme& other) {};
    Elbme& operator=(const Elbme& other) {};
public:
    Elbme(vector<Point*> verts, double L) {
        _verts = verts;
        _L = L;
        _D = _verts.size() - 1;
        _partition.push_back(new Subsimplex(_verts, _L));
        _accuracy = 1e-1;    // should be proportional to problem dimension
        _max_iter = 10;
    };
    vector<Point*> _verts;   // Simplex vertexes used for initialization
    double _L;
    int _D;
    vector<Subsimplex*> _partition;
    vector<Point*> _points;  // Cashe of Point objects to be deleted
    double _accuracy;        // Stopping criteria (tolerance threshold)
    int _max_iter;

    double get_lb_value(vector<Point*> verts, double L, Point* middle_point) {
        double lb_value = -numeric_limits<double>::max();  // Note check numeric_limits<double>min() usage, because the value is 2e-308 ~ 0, and not negative.
        double dist;
        int min_i = numeric_limits<int>::min();
        // 2e-308 value is calculated here, is it correct?
        for (int i=0; i < verts.size(); i++) {
            dist = l2norm(middle_point, verts[i]);
            // cout << lb_value << " is greater thant " << verts[i]->_values[0] - L*dist << endl;
            if (lb_value < verts[i]->_values[0] - L*dist) {
                lb_value = verts[i]->_values[0] - L*dist;
                min_i = i;
            };
        };
        // if (lb_value < 0.0000001) {
        //     cout << "Lb value: "<< lb_value << endl;
        //     cout << " calculated:  L " << L << " dist " << dist << endl; //l2norm(middle_point, verts[min_i]) << endl;
        //     cout << " vert value " << endl; // verts[min_i]->_values[0] << endl;
        //     cout << " min i" << min_i << endl; // verts[min_i]->_values[0] << endl;
        //     exit(0);
        // };
        return lb_value;
    };

    vector<Subsimplex*> divide_subsimplex(Subsimplex* subsimplex, string strategy="longest_half") {
        vector<Subsimplex*> divided_subsimplexes;
        if (strategy== "longest_half") {
            // Find middle point
            double c[_D];
            for (int i=0; i < _D; i++) {
                c[i] = (subsimplex->_le_v1->_X[i] + subsimplex->_le_v2->_X[i]) / 2.;
            };

            Point* middle_point = new Point(c, _D);
            _points.push_back(middle_point);
            middle_point->add_value(get_lb_value(_verts, _L, middle_point));

            // Construct two new simplexes using this middle point.
            Subsimplex* left_subsimplex = new Subsimplex(subsimplex->_L);
            Subsimplex* right_subsimplex = new Subsimplex(subsimplex->_L);

            for (int i=0; i < _D + 1; i++){
                // Point* point = _func->get(new Point(triangle[i], n));
                if (subsimplex->_verts[i] != subsimplex->_le_v1){
                    right_subsimplex->add_vertex(subsimplex->_verts[i]);
                } else {
                    right_subsimplex->add_vertex(middle_point);
                };
                if (subsimplex->_verts[i] != subsimplex->_le_v2) {
                    left_subsimplex->add_vertex(subsimplex->_verts[i]);
                } else {
                    left_subsimplex->add_vertex(middle_point);
                };
            };
            left_subsimplex->init_parameters();
            right_subsimplex->init_parameters();

            divided_subsimplexes.push_back(left_subsimplex);
            divided_subsimplexes.push_back(right_subsimplex);
        };
        return divided_subsimplexes;
    };

    // int get_simplex_to_divide_id(vector<Subsimplex*> subsimplexes) {
    //     // Selects subsimplex with lowest estimated lower bound minimum
    //     int simplex_with_min_lb_value_id = 0;
    //     for (int i=0; i < subsimplexes.size(); i++) {
    //         if (subsimplexes[simplex_with_min_lb_value_id]->_min_lb_value < subsimplexes[i]->_min_lb_value) {
    //             simplex_with_min_lb_value_id = i;
    //         };
    //     };
    //     return simplex_with_min_lb_value_id;
    // };

    Point* minimize() {     // Returns estimated simplex_lb_min_value
        double tolerance = _partition[0]->_tolerance;

        int iter = 0;
        // while (tolerance > _accuracy) {
        while (_max_iter > iter) {
            // Pop subsimplex whith lowest min lb value from partition
            // int simplex_to_divide_id = get_simplex_to_divide_id(_partition);
            Subsimplex* to_divide = _partition[0];
            _partition.erase(_partition.begin() + 0);

            // Divide it
            vector<Subsimplex*> divided_subsimplexes = divide_subsimplex(to_divide);
            // Erase to_divide
            delete to_divide;

            // insert new simplexes as bubble sort
            // Iterate through
            for (int i=0; i < divided_subsimplexes.size(); i++) { // Note: structure which has instert
                _partition.push_back(divided_subsimplexes[i]);    // command would allow to insert as bubble sort
            };

            // sort partition by min lb value
            sort(_partition.begin(), _partition.end(), Subsimplex::compare_by_min_lb_value);

            // cout << "Sorted values: ";
            // for (int i=0; i < _partition.size(); i++) {
            //     cout << _partition[i]->_min_lb_value << ", ";
            // };
            // cout << endl;


            tolerance = _partition[0]->_tolerance;
            iter += 1;
            // cout << it << ". tol " << tolerance << " acc " << _accuracy << " estimate " << _partition[0]->_min_vert_value << " min_lb " << _partition[0]->_min_lb_value << " diameter " << _partition[0]->_diameter << endl;
        };
        return _partition[0]->_min_vert;
    };

    virtual ~Elbme() {
        for (int i=0; i < _points.size(); i++) {
            delete _points[i];
        };
        _points.clear();
        for (int i=0; i < _partition.size(); i++) {
            delete _partition[i];
        };
        _partition.clear();
    };
};

#endif