add an example of 9000 clauses for 3 sat progra
[charm.git] / examples / charm++ / satisfiability / Solver.C
1 #include "main.decl.h"
2 #include "SolverTypes.h"
3 #include "Solver.h"
4 #include <map>
5 #include <vector>
6
7 using namespace std;
8
9 extern CProxy_Main mainProxy;
10 extern int grainsize;
11
12 SolverState::SolverState()
13 {
14     unsolved_clauses = 0; 
15     var_frequency = 0;
16 }
17
18 inline int SolverState::nVars()
19 {
20     return var_frequency; 
21 }
22
23 int SolverState::clausesSize()
24 {
25     return clauses.size();
26 }
27 Var SolverState::newVar(bool sign, bool dvar)
28 {
29     int v = nVars();
30     var_frequency++;
31     return v;
32 }
33
34 void SolverState::attachClause(Clause& c)
35 {
36
37 }
38
39
40 int SolverState::unsolvedClauses()
41 {
42     int unsolved = 0;
43     for(int i=0; i< clauses.size(); i++)
44     {
45         if(clauses[i].size() > 0)
46             unsolved++;
47     }
48
49     return unsolved;
50 }
51
52 void SolverState::printSolution()
53 {
54     for(int _i=0; _i<var_size; _i++)
55     {
56         if(occurrence[_i] == -2)
57             CkPrintf(" TRUE ");
58         else if(occurrence[_i] == -1)
59             CkPrintf(" FALSE ");
60         else
61             CkPrintf(" Anything ");
62     }
63  
64 }
65
66 /* add clause, before adding, check unit conflict */
67 bool SolverState::addClause(CkVec<Lit>& ps)
68 {
69     /*TODO precheck is needed here */
70     if(ps.size() == 1)//unit clause
71     {
72         for(int _i=0; _i<unit_clause_index.size(); _i++)
73         {
74             int index = unit_clause_index[_i];
75
76             Clause unit = clauses[index];
77             Lit     unit_lit = unit[0];
78
79             if(unit_lit == ~ps[0])
80             {
81    
82                 CkPrintf("clause conflict between %d and %d\n", index, clauses.size());
83                 return false;
84             }
85             /*else if(unit_lit == ps[0])
86             {
87                 return true;
88             }*/
89         }
90        /* all unit clauses are checked, no repeat, no conflict */
91         unit_clause_index.push_back(clauses.size());
92     }else{
93     //check whether the clause is already satisfiable in any case
94     /* if x and ~x exist in the same clause, this clause is already satisfiable, we do not have to deal with it */
95    
96     for(int i=0; i< ps.size(); i++)
97     {
98         Lit lit = ps[i];
99
100         for(int j=0; j<i; j++)
101         {
102             if(lit == ~ps[j])
103             {
104                 CkPrintf("This clause is already satisfiable\n");
105                 for(int _i=0; _i<ps.size(); _i++)
106                 {
107                    occurrence[abs(toInt(ps[_i])-1)]--; 
108                    if(toInt(ps[_i]) > 0)
109                        positive_occurrence[abs(toInt(ps[_i])-1)]--; 
110                 }
111                 return true;
112             }
113         }
114     }
115     }
116     clauses.push_back(Clause());
117     clauses[clauses.size()-1].attachdata(ps, false);
118     unsolved_clauses++;
119     // build the linklist for each literal pointing to the clause, where the literal occurs
120     for(int i=0; i< ps.size(); i++)
121     {
122         Lit lit = ps[i];
123
124         if(toInt(lit) > 0)
125             whichClauses[2*toInt(lit)-2].push_back(clauses.size()-1);
126         else
127             whichClauses[-2*toInt(lit)-1].push_back(clauses.size()-1);
128
129     }
130
131     return true;
132 }
133
134 void SolverState::assignclause(CkVec<Clause>& cls )
135 {
136     clauses.removeAll();
137     for(int i=0; i<cls.size(); i++)
138     {
139         clauses.push_back( cls[i]);
140     }
141 }
142
143 /* *********** Solver chare */
144
145 Solver::Solver(SolverState* state_msg)
146 {
147
148     /* Which variable get assigned  */
149     Lit lit = state_msg->assigned_lit;
150 #ifdef DEBUG    
151     CkPrintf("\n\nNew chare: literal = %d, occurrence size=%d, level=%d \n", toInt(lit), state_msg->occurrence.size(), state_msg->level);
152 #endif    
153     SolverState *next_state = copy_solverstate(state_msg);
154     
155     //Unit clauses
156     /* use this value to propagate the clauses */
157     // deal with the clauses where this literal is
158     int _unit_ = -1;
159    
160 unit_propagation:while(1){
161     int pp_ = 1;
162     int pp_i_ = 2;
163     int pp_j_ = 1;
164
165     if(toInt(lit) < 0)
166     {
167         pp_ = -1;
168         pp_i_ = 1;
169         pp_j_ = 2;
170     }
171
172     next_state->occurrence[pp_*toInt(lit)-1] = -pp_i_;
173     
174     //CkPrintf(" index=%d, total size=%d, original msg size=%d\n", pp_*2*toInt(lit)-pp_i_, next_state->whichClauses.size(), state_msg->whichClauses.size());
175
176     CkVec<int> &inClauses = next_state->whichClauses[pp_*2*toInt(lit)-pp_i_];
177     CkVec<int> &inClauses_opposite = next_state->whichClauses[pp_*2*toInt(lit)-pp_j_];
178
179     // literal with same sign, remove all these clauses
180     for(int j=0; j<inClauses.size(); j++)
181     {
182         int cl_index = inClauses[j];
183
184         Clause& cl_ = next_state->clauses[cl_index];
185         //for all the literals in this clauses, the occurrence decreases by 1
186         for(int k=0; k< cl_.size(); k++)
187         {
188             Lit lit_ = cl_[k];
189             if(toInt(lit_) == toInt(lit))
190                 continue;
191             next_state->occurrence[abs(toInt(lit_)) - 1]--;
192             if(toInt(lit_) > 0)
193             {
194                 next_state->positive_occurrence[toInt(lit_)-1]--;
195                 //remove the clause index for the literal
196                 for(int _i = 0; _i<next_state->whichClauses[2*toInt(lit_)-2].size(); _i++)
197                 {
198                     if(next_state->whichClauses[2*toInt(lit_)-2][_i] == cl_index)
199                     {
200                         next_state->whichClauses[2*toInt(lit_)-2].remove(_i);
201                         break;
202                     }
203                 }
204             }else
205             {
206                 for(int _i = 0; _i<next_state->whichClauses[-2*toInt(lit_)-1].size(); _i++)
207                 {
208                     if(next_state->whichClauses[-2*toInt(lit_)-1][_i] == cl_index)
209                     {
210                         next_state->whichClauses[-2*toInt(lit_)-1].remove(_i);
211                         break;
212                     }
213                 }
214             }
215         }
216             next_state->clauses[cl_index].resize(0);
217
218     }
219    
220     for(int j=0; j<inClauses_opposite.size(); j++)
221     {
222         int cl_index_ = inClauses_opposite[j];
223         Clause& cl_neg = next_state->clauses[cl_index_];
224         cl_neg.remove(-toInt(lit));
225         
226         //becomes a unit clause
227          if(cl_neg.size() == 1)
228          {
229              next_state->unit_clause_index.push_back(cl_index_);
230          }else if (cl_neg.size() == 0)
231          {
232                 return;
233          }
234     }
235
236
237     _unit_++;
238     if(_unit_ == next_state->unit_clause_index.size())
239         break;
240     Clause cl = next_state->clauses[next_state->unit_clause_index[_unit_]];
241
242     while(cl.size() == 0){
243         _unit_++;
244         if(_unit_ == next_state->unit_clause_index.size())
245             break;
246         cl = next_state->clauses[next_state->unit_clause_index[_unit_]];
247         
248     };
249
250     if(_unit_ == next_state->unit_clause_index.size())
251         break;
252     lit = cl[0];
253     }
254     /***************/
255
256     int unsolved = next_state->unsolvedClauses();
257     if(unsolved == 0)
258     {
259         CkPrintf("One solution found in parallel processing \n");
260         
261 #ifdef SOLUTION
262         next_state->printSolution();
263 #endif
264
265         mainProxy.done();
266         return;
267     }
268     int max_index = get_max_element(next_state->occurrence);
269 #ifdef DEBUG
270     CkPrintf("max index = %d\n", max_index);
271 #endif
272     //if() left literal unassigned is larger than a grain size, parallel 
273     ///* When we start sequential 3SAT Grainsize problem*/
274    
275     /* the other branch */
276     SolverState *new_msg2 = copy_solverstate(next_state);;
277     new_msg2->var_size = state_msg->var_size;
278
279     next_state->level = state_msg->level+1;
280
281     int lower = state_msg->lower;
282     int higher = state_msg->higher;
283     int middle = (lower+higher)/2;
284     int positive_max = next_state->positive_occurrence[max_index];
285     if(positive_max >= next_state->occurrence[max_index] - positive_max)
286     {
287         next_state->assigned_lit = Lit(max_index+1);
288         next_state->occurrence[max_index] = -2;
289     }
290     else
291     {
292         next_state->assigned_lit = Lit(-max_index-1);
293         next_state->occurrence[max_index] = -1;
294     }
295     bool satisfiable_1 = true;
296     bool satisfiable_0 = true;
297
298     if(unsolved > grainsize)
299     {
300         next_state->lower = lower + 1;
301         next_state->higher = middle;
302         *((int *)CkPriorityPtr(next_state)) = lower+1;
303         CkSetQueueing(next_state, CK_QUEUEING_IFIFO);
304         CProxy_Solver::ckNew(next_state);
305     }
306     else //sequential
307     {
308         satisfiable_1 = seq_solve(next_state);
309         if(satisfiable_1)
310         {
311             CkPrintf(" One solutions found by sequential algorithm\n");
312             mainProxy.done();
313             return;
314         }
315     }
316
317     new_msg2->level = state_msg->level+1;
318     if(positive_max >= new_msg2->occurrence[max_index] - positive_max)
319     {
320         new_msg2->assigned_lit = Lit(-max_index-1);
321         new_msg2->occurrence[max_index] = -1;
322     }
323     else
324     {
325         new_msg2->assigned_lit = Lit(max_index+1);
326         new_msg2->occurrence[max_index] = -2;
327     }
328     unsolved = new_msg2->unsolvedClauses();
329     if(unsolved > grainsize)
330     {
331
332         new_msg2->lower = middle + 1;
333         new_msg2->higher = higher-1;
334         *((int *)CkPriorityPtr(new_msg2)) = middle+1;
335         CkSetQueueing(new_msg2, CK_QUEUEING_IFIFO);
336         CProxy_Solver::ckNew(new_msg2);
337     }
338     else
339     {
340         satisfiable_0 = seq_solve(new_msg2);
341
342         if(satisfiable_0)
343         {
344             CkPrintf(" One solutions found by sequential algorithm\n");
345             mainProxy.done();
346             return;
347         }
348
349     }
350
351 }
352
353 /* Which literals are already assigned, which is assigned this interation, the unsolved clauses */
354 /* should all these be passed as function parameters */
355 /* solve the 3sat in sequence */
356
357 long long int computes = 0;
358 bool Solver::seq_solve(SolverState* state_msg)
359 {
360     /* Which variable get assigned  */
361     Lit lit = state_msg->assigned_lit;
362        
363 #ifdef DEBUG
364     CkPrintf("\n\n Computes=%d Sequential SAT New chare: literal = %d,  level=%d, unsolved clauses=%d\n", computes++, toInt(assigned_var), state_msg->level, state_msg->clauses.size());
365     //CkPrintf("\n\n Computes=%d Sequential SAT New chare: literal = %d, occurrence size=%d, level=%d \n", computes++, toInt(assigned_var), state_msg->occurrence.size(), state_msg->level);
366 #endif
367     SolverState *next_state = copy_solverstate(state_msg);
368     
369     //Unit clauses
370     /* use this value to propagate the clauses */
371 #ifdef DEBUG
372     CkPrintf(" remainning clause size is %d\n", (state_msg->clauses).size());
373 #endif
374
375     int _unit_ = -1;
376     while(1){
377     int pp_ = 1;
378     int pp_i_ = 2;
379     int pp_j_ = 1;
380
381     if(toInt(lit) < 0)
382     {
383         pp_ = -1;
384         pp_i_ = 1;
385         pp_j_ = 2;
386     }
387
388     next_state->occurrence[pp_*toInt(lit)-1] = -pp_i_;
389     
390     CkVec<int> &inClauses = next_state->whichClauses[pp_*2*toInt(lit)-pp_i_];
391     CkVec<int> &inClauses_opposite = next_state->whichClauses[pp_*2*toInt(lit)-pp_j_];
392
393     // literal with same sign, remove all these clauses
394     for(int j=0; j<inClauses.size(); j++)
395     {
396         int cl_index = inClauses[j];
397
398         Clause& cl_ = next_state->clauses[cl_index];
399         //for all the literals in this clauses, the occurrence decreases by 1
400         for(int k=0; k< cl_.size(); k++)
401         {
402             Lit lit_ = cl_[k];
403             if(toInt(lit_) == toInt(lit))
404                 continue;
405             next_state->occurrence[abs(toInt(lit_)) - 1]--;
406             if(toInt(lit_) > 0)
407             {
408                 next_state->positive_occurrence[toInt(lit_)-1]--;
409                 //remove the clause index for the literal
410                 for(int _i = 0; _i<next_state->whichClauses[2*toInt(lit_)-2].size(); _i++)
411                 {
412                     if(next_state->whichClauses[2*toInt(lit_)-2][_i] == cl_index)
413                     {
414                         next_state->whichClauses[2*toInt(lit_)-2].remove(_i);
415                         break;
416                     }
417                 }
418             }else
419             {
420                 for(int _i = 0; _i<next_state->whichClauses[-2*toInt(lit_)-1].size(); _i++)
421                 {
422                     if(next_state->whichClauses[-2*toInt(lit_)-1][_i] == cl_index)
423                     {
424                         next_state->whichClauses[-2*toInt(lit_)-1].remove(_i);
425                         break;
426                     }
427                 }
428             }
429         }
430             next_state->clauses[cl_index].resize(0);
431     }
432    
433     for(int j=0; j<inClauses_opposite.size(); j++)
434     {
435         int cl_index_ = inClauses_opposite[j];
436         Clause& cl_neg = next_state->clauses[cl_index_];
437         cl_neg.remove(-toInt(lit));
438             //becomes a unit clause
439          if(cl_neg.size() == 1)
440          {
441              next_state->unit_clause_index.push_back(cl_index_);
442          }else if (cl_neg.size() == 0)
443          {
444     //            CkPrintf(" conflict found in seqential testing!\n");
445                 return false;
446          }
447     }
448    
449     _unit_++;
450     if(_unit_ == next_state->unit_clause_index.size())
451         break;
452     Clause cl = next_state->clauses[next_state->unit_clause_index[_unit_]];
453     
454     while(cl.size() == 0){
455         _unit_++;
456         if(_unit_ == next_state->unit_clause_index.size())
457             break;
458         cl = next_state->clauses[next_state->unit_clause_index[_unit_]];
459         
460     };
461
462     if(_unit_ == next_state->unit_clause_index.size())
463         break;
464     
465     lit = cl[0];
466     }
467    
468    
469     int unsolved = next_state->unsolvedClauses();
470     if(unsolved == 0)
471     {
472         CkPrintf("One solution found\n");
473 #ifdef SOLUTION
474         next_state->printSolution();
475 #endif
476         return true;
477     }
478     
479     /**********************/
480     
481         /* it would be better to insert the unit literal in order of their occurrence */
482         /* make a decision and then fire new tasks */
483         /* if there is unit clause, should choose these first??? TODO */
484         /* TODO which variable to pick up */
485         /*unit clause literal and also which occurrs most times */
486         int max_index =  get_max_element(next_state->occurrence);
487 #ifdef DEBUG
488         CkPrintf("max index = %d\n", max_index);
489 #endif
490         next_state->var_size = state_msg->var_size;
491         next_state->level = state_msg->level+1;
492
493         int positive_max = next_state->positive_occurrence[max_index];
494         if(positive_max >= next_state->occurrence[max_index] - positive_max)
495         {
496             next_state->occurrence[max_index] = -2;
497             next_state->assigned_lit = Lit(max_index+1);
498         }
499         else
500         {
501             next_state->occurrence[max_index] = -1;
502             next_state->assigned_lit = Lit(-max_index-1);
503         } 
504
505         bool   satisfiable_1 = seq_solve(next_state);
506         if(satisfiable_1)
507         {
508             return true;
509         }
510         
511        
512         if(positive_max >= next_state->occurrence[max_index] - positive_max)
513         {
514             next_state->occurrence[max_index] = -1;
515             next_state->assigned_lit = Lit(-max_index-1);
516         }
517         else
518         {
519             next_state->occurrence[max_index] = -2;
520             next_state->assigned_lit = Lit(max_index+1);
521         } 
522             
523         bool satisfiable_0 = seq_solve(next_state);
524         if(satisfiable_0)
525         {
526             return true;
527         }
528
529         //CkPrintf("Unsatisfiable through sequential\n");
530         return false;
531
532 }