Project

General

Profile

Statistics
| Branch: | Revision:

cool / src / lib / GMLMIP-0.1 / formulas / formula.c @ 7c4d2eb4

History | View | Annotate | Download (10.2 KB)

1
#include "formula.h"
2

    
3
#include <set>
4
#include <cassert>
5

    
6
IFormula::~IFormula() {
7
}
8

    
9
template<class ModalValueType>
10
Formula<ModalValueType>::Formula(){
11
        number_of_variables = 0;
12
        variables.set_empty_key(-1);            //
13
        variables_back.set_empty_key(-1);       //
14
        data d;
15
        modal_variables.set_empty_key(d); 
16
        modal_variables_back.set_empty_key(-1); //
17
        Premise<ModalValueType> p;
18
        rule_cache.set_empty_key(p);
19
}
20

    
21
template<class ModalValueType>
22
bdd Formula<ModalValueType>::variable(int n){
23
        if(variables.count(n) == 0){
24
                variables[n] = number_of_variables;
25
                variables_back[number_of_variables] = n;
26
                number_of_variables++;
27
        }
28
        return bdd_ithvarpp(variables[n]);
29
}
30

    
31
template<class ModalValueType>
32
bdd Formula<ModalValueType>::modal_var(bdd *b, ModalValueType &n){
33
        data d(n, b);
34
        if(modal_variables.count(d) == 0){
35
                modal_variables[d] = number_of_variables;
36
                modal_variables_back[number_of_variables] = d;
37
                number_of_variables++;
38
        }
39
        return bdd_ithvarpp(modal_variables[d]);
40
}
41

    
42

    
43
template<class ModalValueType>
44
void Formula<ModalValueType>::print_back_vars(){
45
        cout << setw(10) << "BDD number   " << setw(10) << "Input Number" << endl;
46
        for(int i=0; i < 1000; i++)
47
                if(variables_back.count(i) > 0 )
48
                        cout <<  setw(10) << i << setw(10) << variables_back[i] << endl;
49
        cout << endl;
50
}
51

    
52
vector<SatisfyingAssignment> gatherer_stack;
53

    
54
void sat_handler_gatherer(char *varset, int size){
55
        SatisfyingAssignment s(size, varset);
56
        gatherer_stack.push_back(s);
57
}
58

    
59
template<class ModalValueType>
60
set<set<int> > Formula<ModalValueType>::gatherRules(bdd b) {
61
        bdd_allsat(b, sat_handler_gatherer);
62
        vector<SatisfyingAssignment> sat_ass = gatherer_stack;
63
        gatherer_stack.clear();
64
        set<set<int> > disj;
65
        for (vector<SatisfyingAssignment>::iterator it = sat_ass.begin();
66
                it != sat_ass.end(); ++it)
67
        {
68
                set<int> conj;
69
                SatisfyingAssignment& sat = *it;
70
                for (int v = 0; v < sat.get_size(); v++) {
71
                        if (sat.get_array_i(v) != 1 && sat.get_array_i(v) != 0) continue;
72
                        if (!variables_back.count(v)) continue;
73
                        int sgn = ((sat.get_array_i(v)==1) ? 1 : (-1));
74
                        conj.insert(variables_back[v] * sgn);
75
                }
76
                disj.insert(conj);
77
        }
78
        return disj;
79
}
80

    
81
void sat_handler(char *varset, int size){
82
        SatisfyingAssignment s(size, varset);
83
        sat_ass_stack.push_back(s);
84
}
85

    
86
template<class ModalValueType>
87
RuleCollection Formula<ModalValueType>::check_satisfiability(bdd b){
88
        //bdd_printdot(b);
89
        bdd_allsat(b, sat_handler); // pushes satisfying assignment onto stack
90
        vector<SatisfyingAssignment> sat_ass = sat_ass_stack; 
91
        sat_ass_stack.clear();        
92
        // it might be zero if already on the syntactic level there is something like A /\ ~A
93
        assert(sat_ass.size() <= 1);
94
        if (sat_ass.size() == 1) {
95
                return apply_rules(sat_ass[0]);
96
        } else {
97
                // TODO: return {} or { {} }?
98
                // the secound one would mean { False }
99
                return set<set<set<int> > >();
100
        }
101
}
102

    
103
template<class ModalValueType>
104
RuleCollection Formula<ModalValueType>::apply_rules(SatisfyingAssignment& sat){
105
        RuleCollection result;
106
        
107
        int no_of_positive = 0;
108
        int no_of_negative = 0;
109
        bdd **underlying_formulae;
110
        ModalValueType *positive_modal_indices, *negative_modal_indices;
111

    
112
        // We first find the size of the required arrays.
113
        for(int v=0; v < sat.get_size(); v++){
114
                if(sat.get_array_i(v)==1 && modal_variables_back.count(v)) // if true and represents a modal formula.
115
                        no_of_positive++;
116
                else if(sat.get_array_i(v)==0 && modal_variables_back.count(v)) // if false and ...
117
                        no_of_negative++;
118
        }
119
        // cout << no_of_positive + no_of_negative << endl;
120
        if(no_of_positive + no_of_negative > 0){
121
        
122
                // We then initialize our new arrays.
123
                underlying_formulae = new bdd* [no_of_positive+no_of_negative];
124
        
125
                positive_modal_indices = new ModalValueType[no_of_positive];
126
                negative_modal_indices = new ModalValueType[no_of_negative];
127

    
128
                // We then fill our underlying formulae and modal value arrays.
129
                int i=0, j=0;
130
                for(int v=0; v < sat.get_size(); v++){
131
                        if(sat.get_array_i(v)==1 && modal_variables_back.count(v)){
132
                                underlying_formulae[i] = (modal_variables_back[v]).argument;
133
                                positive_modal_indices[i] = (modal_variables_back[v]).value;
134
                                i++;
135
                        }
136
                        else if(sat.get_array_i(v)==0 && modal_variables_back.count(v)){ //if false and represents a modal formula.
137
                                underlying_formulae[no_of_positive+j] = (modal_variables_back[v]).argument;
138
                                negative_modal_indices[j] = (modal_variables_back[v]).value;
139
                                j++;
140
                        }
141
                }
142

    
143
                // Construct a GML premise
144
                Premise<ModalValueType> premise(no_of_positive, no_of_negative, positive_modal_indices, negative_modal_indices);
145
                // Check the rule cache
146
                if(rule_cache.count(premise)){ // if in cache
147
                        SetOfConclusions concs = rule_cache[premise];
148
                        for(int i = 0; i < concs.number_of_conclusions(); i++){
149
                                bdd b = concs.get_jth_conclusion(underlying_formulae, i);
150
                                result.insert(gatherRules(b));
151
                        }
152
                } else {
153
                        RuleCollection r = enumerate_rules(premise, underlying_formulae);
154
                        result.insert(r.begin(), r.end());
155
                }
156

    
157
                delete [] underlying_formulae;
158
                delete [] positive_modal_indices;
159
                delete [] negative_modal_indices; 
160
        }
161
        
162
        return result;
163
}
164

    
165
template<class ModalValueType>
166
RuleCollection Formula<ModalValueType>::enumerate_rules(Premise<ModalValueType>& prem, bdd** underlying_formulae){
167
        RuleCollection result;
168
        vector<vector<bool> > valid_nodes;
169
        vector<pair<vector<bool>, int> > anti_nodes;
170
        vector<vector<bool> > node_queue;
171
        int total_valuations = prem.get_total_valuations();
172

    
173
        glp_prob *problem = glp_create_prob();
174

    
175
        // Parameters while solving with glpk
176
        glp_iocp* parameters = new glp_iocp;
177
        glp_init_iocp(parameters);
178
        parameters->presolve = GLP_ON;
179

    
180
        parameters->msg_lev = GLP_MSG_ERR;
181
        //parameters->msg_lev = GLP_MSG_ALL;
182
        
183
        load_linear_program(problem, prem);
184
        
185
        vector<bool> zero(total_valuations, false);
186
        
187
        if (run_solver(problem, parameters, false, zero, -1)) {
188
                if (run_solver(problem, parameters, true, zero, 0)) {
189
                } else {
190
                        node_queue.push_back(zero);
191
                        //int i = 0;
192
                        while(node_queue.size()) {
193
                                //cout << "outer space: " << i++ << "  runs: " << test_counter << "  queue: " << node_queue.size() << "  sols: " << valid_nodes.size() << "  antisols: " << anti_nodes.size() << endl;
194
                                RuleCollection r = enum_rules_intern(problem, parameters, prem, underlying_formulae, node_queue, valid_nodes, anti_nodes);
195
                                result.insert(r.begin(), r.end());
196
                        }
197
                        if (!rule_cache.count(prem)) {
198
                                rule_cache[prem] = SetOfConclusions(prem.get_n() + prem.get_m(), valid_nodes);
199
                        }
200
                        //cout << "runs: " << test_counter << "  sols: " << valid_nodes.size() << "  antisols: " << anti_nodes.size() << endl;
201
                }
202
        }
203

    
204
        //Clean-up
205
        delete parameters;
206
        glp_delete_prob(problem);
207
        
208
        return result;
209
}
210

    
211
template<class ModalValueType>
212
RuleCollection Formula<ModalValueType>::enum_rules_intern(glp_prob* problem, glp_iocp* parameters, Premise<ModalValueType>& prem, bdd** underlying_formulae,
213
                                                                                                vector<vector<bool> >& node_queue, vector<vector<bool> >& valid_nodes, vector<pair<vector<bool>, int> >& anti_nodes){
214
        RuleCollection result;
215
        vector<vector<bool> > new_node_queue;
216
        int total_valuations = prem.get_total_valuations();
217

    
218
        for (unsigned int k = 0; k < node_queue.size(); k++) {
219
                vector<bool> current_node = node_queue[k];
220

    
221
                int power = total_valuations-1;
222
                while(power >= 0 && !current_node[power])
223
                        power--;
224
                for(int i = power+1; i < total_valuations; i++){
225
                        vector<bool> new_node(current_node);
226
                        new_node[i] = true;
227

    
228
                        // Check node isn't superset of something already valid
229
                        bool isSubnode = false;
230
                        for(unsigned int i1=0; i1 < valid_nodes.size() && !isSubnode; i1++){
231
                                bool isSub = true;
232
                                for(int j1=0; j1 < total_valuations; j1++){
233
                                        if (valid_nodes[i1][j1] && !new_node[j1]) {
234
                                                isSub = false;
235
                                                break;
236
                                        }
237
                                }
238
                                if (isSub) isSubnode = true;
239
                        }
240
                        // Check node isn't superset of something of an anti-solution
241
                        for(unsigned int i1=0; i1 < anti_nodes.size() && !isSubnode; i1++){
242
                                bool isSub = true;
243
                                int idx = anti_nodes[i1].second;
244
                                vector<bool> antinode = anti_nodes[i1].first;
245
                                for(int j1=0; j1 < total_valuations; j1++){
246
                                        if ((j1 <= idx && antinode[j1] != new_node[j1]) || (j1 > idx && antinode[j1] && !new_node[j1])) {
247
                                                isSub = false;
248
                                                break;
249
                                        }
250
                                }
251
                                if (isSub) isSubnode = true;
252
                        }
253

    
254
                        if (!isSubnode) {
255
                                if (run_solver(problem, parameters, false, new_node, i)) {
256
                                        if (run_solver(problem, parameters, true, new_node, 0)) {
257
                                                valid_nodes.push_back(new_node);
258
                                                //cout << "solution found " << valid_nodes.size() << endl;
259
                                                //for (int l = 0; l < total_valuations; l++)
260
                                                //        cout << "   " << new_node[l] << endl;
261
                                                vector<vector<bool> > singleton;
262
                                                singleton.push_back(new_node);
263
                                                SetOfConclusions temp = SetOfConclusions(prem.get_n() + prem.get_m(), singleton);
264
                                                bdd b = temp.get_jth_conclusion(underlying_formulae, 0);
265
                                                result.insert(gatherRules(b));
266
                                        } else {
267
                                                new_node_queue.push_back(new_node);
268
                                        }
269
                                }
270
                                else {
271
                                        pair<vector<bool>, int> antipair = pair<vector<bool>, int>(new_node, i);
272
                                        anti_nodes.push_back(antipair);
273
                                }
274
                        }
275
                }
276
        }
277
        node_queue = new_node_queue;
278
        return result;
279
}        
280

    
281
template<class ModalValueType>
282
bool Formula<ModalValueType>::run_solver(glp_prob *problem, glp_iocp* parameters, bool strict, vector<bool> current_node, int index){
283
        test_counter++;
284

    
285
        // Modify glp problem to match node
286
        for(unsigned int i = 0; i < current_node.size(); i++){
287
                if(current_node[i])
288
                        glp_set_row_bnds(problem, i+1, GLP_LO, 1.0, 0.0); // > 0
289
                        //+1 since valuations go zero to 2^(n+m)-1, glpk rows go 1 to 2^(n+m).
290
                else if (strict || static_cast<int>(i) <= index)
291
                        glp_set_row_bnds(problem, i+1, GLP_UP, 0.0, 0.0); // leq 0
292
                else
293
                        glp_set_row_bnds(problem, i+1, GLP_FR, 0.0, 0.0); // don't care
294
        }
295

    
296
        if (!strict) {
297
                parameters->tm_lim = 1000;
298
        } else {
299
                parameters->tm_lim = 1000000000;
300
        }
301

    
302
        // Solve the LP
303
        int result = glp_intopt(problem, parameters);
304
                
305
        if(result == 0 || result == GLP_ENOPFS) { //if no errors from solver or a relaxtion feasibility error
306
                if(glp_mip_status(problem)==GLP_OPT || glp_mip_status(problem)==GLP_FEAS){ // if feasible
307
                        return true;
308
                }
309
                else {
310
                        return false;
311
                }
312
        }
313
        else if (!strict && result == GLP_ETMLIM) { // if the solver timed out for non-strict problems
314
                //cout << "timeout" << endl;
315
                return true;
316
        } else { // if error but not a feasibility error
317
                cout << "glpk error: " << result << endl;
318
                exit(1);
319
        }
320
}
321

    
322

    
323
template class Formula<int>;
324
template class Formula<Rational>;