source: vis_dev/cusp-1.1/src/smt/smtMp.c @ 23

Last change on this file since 23 was 12, checked in by cecile, 13 years ago

cusp added

File size: 14.4 KB
Line 
1/**CFile***********************************************************************
2
3  FileName    [smtMp.c]
4
5  PackageName [smt]
6
7  Synopsis    [Routines for smt mp function.]
8
9  Author      [Hyondeuk Kim]
10
11  Copyright   [Copyright (c) 1995-2004, Regents of the University of Colorado
12
13  All rights reserved.
14
15  Redistribution and use in source and binary forms, with or without
16  modification, are permitted provided that the following conditions
17  are met:
18
19  Redistributions of source code must retain the above copyright
20  notice, this list of conditions and the following disclaimer.
21
22  Redistributions in binary form must reproduce the above copyright
23  notice, this list of conditions and the following disclaimer in the
24  documentation and/or other materials provided with the distribution.
25
26  Neither the name of the University of Colorado nor the names of its
27  contributors may be used to endorse or promote products derived from
28  this software without specific prior written permission.
29
30  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
31  "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
32  LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
33  FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
34  COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
35  INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
36  BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
37  LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
38  CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
39  LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
40  ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
41  POSSIBILITY OF SUCH DAMAGE.]
42
43******************************************************************************/
44
45#ifdef HAVE_LIBGMP
46
47#include "smt.h"
48
49smtMp_t *
50smt_mp_init(smtManager_t * sm)
51{
52  smtMp_t * mp = (smtMp_t *) malloc(sizeof(smtMp_t));
53 
54  sm->mp = mp;
55
56  smt_mp_init_pool(mp);
57
58  mpq_init(mp->zero);
59
60  mpq_init(mp->one);
61
62  mpq_init(mp->m_one);
63
64  mpq_set_d(mp->zero, 0);
65
66  mpq_set_d(mp->one, 1);
67
68  mpq_set_d(mp->m_one, -1);
69
70  smt_mp_init_atom_constants(sm);
71
72  smt_mp_init_nvar_values(sm);
73 
74  return mp;
75}
76
77void
78smt_mp_init_pool(smtMp_t * mp)
79{
80  mpq_t * pool;
81  int i;
82
83  mp->plimit = 10;
84
85  pool = (mpq_t *) malloc(sizeof(mpq_t) * mp->plimit);
86
87  mp->pool = pool;
88
89  for(i = 0; i < mp->plimit; i++) {
90    mpq_init(pool[i]);
91  }
92}
93
94void
95smt_mp_init_nvar_values(smtManager_t * sm)
96{
97  smtMp_t * mp = sm->mp;
98  int i, size;
99
100  mp->num_nvar = sm->nvarArr->size;
101
102  size = mp->num_nvar + 1;
103
104  mp->values = (mpq_t *) malloc(sizeof(mpq_t) * size);
105  mp->pvalues = (mpq_t *) malloc(sizeof(mpq_t) * size);
106
107  for(i = 0; i < size; i++) {
108    mpq_init(mp->values[i]);
109    mpq_init(mp->pvalues[i]);
110    mpq_set_d(mp->values[i], 0);
111    mpq_set_d(mp->pvalues[i], 0);
112  }
113}
114
115void
116smt_mp_init_atom_constants(smtManager_t * sm)
117{
118  smtMp_t * mp = sm->mp;
119  smtAvar_t * avar;
120  int i, size, id;
121
122  mp->num_avar = sm->avarArr->size;
123 
124  size = mp->num_avar + 1;
125
126  mp->constants = (mpq_t *) malloc(sizeof(mpq_t) * size);
127  memset(mp->constants, 0, sizeof(mpq_t) * size);
128  mp->weights = (mpq_t *) malloc(sizeof(mpq_t) * size);
129  memset(mp->weights, 0, sizeof(mpq_t) * size);
130
131  for(i = 0; i < sm->avarArr->size; i++) {
132    avar = (smtAvar_t *) sm->avarArr->space[i];
133    id = avar->id;
134    mpq_init(mp->constants[id]);
135    mpq_init(mp->weights[id]);
136    mpq_set_d(mp->weights[id], 0);
137  }
138
139  return;
140}
141
142void
143smt_mp_free(smtMp_t * mp)
144{
145  smt_mp_free_pool(mp);
146
147  mpq_clear(mp->zero);
148
149  mpq_clear(mp->one);
150
151  mpq_clear(mp->m_one);
152
153  smt_mp_free_atom_constants(mp);
154
155  smt_mp_free_nvar_values(mp);
156 
157  free(mp);
158}
159
160void
161smt_mp_free_pool(smtMp_t * mp)
162{
163  mpq_t * pool = mp->pool;
164  int i;
165
166  for(i = 0; i < mp->plimit; i++) {
167    mpq_clear(pool[i]);
168  }
169
170  free(pool);
171}
172
173void
174smt_mp_free_atom_constants(smtMp_t * mp)
175{
176  int i, id;
177
178  for(i = 0; i < mp->num_avar; i++) {
179    id = i + 1;
180    mpq_clear(mp->constants[id]);
181    mpq_clear(mp->weights[id]);
182  }
183
184  free(mp->constants);
185  free(mp->weights);
186
187  return;
188}
189
190void
191smt_mp_free_nvar_values(smtMp_t * mp)
192{
193  int i, size;
194
195  size = mp->num_nvar + 1;
196
197  for(i = 0; i < size; i++) {
198    mpq_clear(mp->values[i]);
199    mpq_clear(mp->pvalues[i]);
200  }
201
202  free(mp->values);
203  free(mp->pvalues);
204}
205
206void
207smt_mp_assign_atom_constants(smtManager_t * sm)
208{
209  smtMp_t * mp = sm->mp;
210  smtFml_t * avfml;
211  smtAvar_t * avar;
212  int i;
213
214  for(i = 0; i < sm->avfmlArr->size; i++) {
215    avfml = (smtFml_t *) sm->avfmlArr->space[i];
216    avar = avfml->avar;
217    if (!avar) continue;
218    smt_mp_assign_atom_constant(mp, avfml);   
219  }
220
221  return;
222}
223
224void
225smt_mp_assign_atom_constant(smtMp_t * mp, smtFml_t * avfml)
226{
227  mpq_t * pool = mp->pool;
228  mpq_t * mp_constant = 0, * coeff = 0;
229  smtAvar_t * avar = avfml->avar;
230  smtFml_t * lfml, * rfml;
231  int nminus;
232
233  /* mp vars */
234  mp_constant = &pool[0];
235  coeff = &pool[1];
236 
237  lfml = (smtFml_t *) avfml->subfmlArr->space[0];
238  rfml = (smtFml_t *) avfml->subfmlArr->space[1];
239   
240  mpq_set_d(*mp_constant, 0);
241 
242  mpq_set_d(*coeff, 1);
243  nminus = 0;
244 
245  smt_mp_get_atom_constant(mp, lfml, mp_constant, coeff, nminus);
246 
247  mpq_set_d(*coeff, 1);
248  nminus = 1;
249 
250  smt_mp_get_atom_constant(mp, rfml, mp_constant, coeff, nminus);
251
252  /* assign atom constant */
253  mpq_set(mp->constants[avar->id], *mp_constant);
254 
255  return;
256}
257
258void
259smt_mp_get_atom_constant(
260  smtMp_t * mp,
261  smtFml_t * fml,
262  mpq_t * mp_constant,
263  mpq_t * coeff,
264  int nminus)
265{
266  mpq_t * pool = mp->pool;
267  mpq_t * tmp_a = 0, * tmp_b = 0;
268  smtFml_t * lfml, * rfml, * tfml;
269  char * str_a, * str_b;
270  int i;
271
272  /* mp vars */
273  tmp_a = &pool[2];
274  tmp_b = &pool[3];
275 
276  if (fml->type == NUM_c) {
277    mpq_set_str(*tmp_a, (char *) fml->subfmlArr->space[0], 10);
278    if (nminus%2 == 1) {
279      mpq_mul(*tmp_a, *coeff, *tmp_a);
280    } else {
281      mpq_neg(*tmp_a, *tmp_a);
282      mpq_mul(*tmp_a, *coeff, *tmp_a);
283    }
284
285    mpq_add(*mp_constant, *mp_constant, *tmp_a);
286   
287  } else if (fml->type == MINUS_c) {
288    nminus++;
289    tfml = (smtFml_t *) fml->subfmlArr->space[0];
290    smt_mp_get_atom_constant(mp, tfml, mp_constant, coeff, nminus);
291   
292  } else if (fml->type == SUB_c) {
293    lfml = (smtFml_t *) fml->subfmlArr->space[0];
294    rfml = (smtFml_t *) fml->subfmlArr->space[1];
295   
296    smt_mp_get_atom_constant(mp, lfml, mp_constant, coeff, nminus);
297    smt_mp_get_atom_constant(mp, rfml, mp_constant, coeff, nminus + 1);
298
299  } else if (fml->type == ADD_c) {
300    for(i = 0; i < fml->subfmlArr->size; i++) {
301      tfml = (smtFml_t *) fml->subfmlArr->space[i];
302      smt_mp_get_atom_constant(mp, tfml, mp_constant, coeff, nminus);
303    }
304   
305  } else if (fml->type == MUL_c) {
306    lfml = (smtFml_t *) fml->subfmlArr->space[0];
307    rfml = (smtFml_t *) fml->subfmlArr->space[1];   
308   
309    if (lfml->type == NUM_c && rfml->type == NUM_c) {
310      mpq_set_str(*tmp_a, (char *) lfml->subfmlArr->space[0], 10);
311      mpq_set_str(*tmp_b, (char *) rfml->subfmlArr->space[0], 10);
312
313      mpq_mul(*tmp_b, *tmp_a, *tmp_b);
314       
315      if(nminus%2 == 1) {
316        mpq_mul(*tmp_b, *tmp_b, *coeff);
317      } else {
318        mpq_neg(*tmp_a, *coeff);
319        mpq_mul(*tmp_b, *tmp_a, *tmp_b);
320      }
321
322      mpq_add(*mp_constant, *mp_constant, *tmp_b);     
323    } else {
324      fprintf(stdout, "ERROR: WRONG LINEAR ARITHMETIC ATOM TYPE\n");
325      str_a = smt_convert_fml_to_string(lfml);
326      str_b = smt_convert_fml_to_string(rfml);
327      fprintf(stdout, "fml = %s * %s\n", str_a, str_b);
328      fprintf(stdout, "unknown\n");
329      exit(0);     
330    }
331
332  } else if (fml->type == DIV_c) {
333    lfml = (smtFml_t *) fml->subfmlArr->space[0];
334    rfml = (smtFml_t *) fml->subfmlArr->space[1];
335   
336    if (lfml->type == NUM_c && rfml->type == NUM_c) {
337      mpq_set_str(*tmp_a, (char *) lfml->subfmlArr->space[0], 10);
338      mpq_set_str(*tmp_b, (char *) rfml->subfmlArr->space[0], 10);
339     
340      mpq_div(*tmp_b, *tmp_a, *tmp_b);
341       
342      if(nminus%2 == 1) {
343        mpq_mul(*tmp_b, *tmp_b, *coeff);
344      } else {
345        mpq_neg(*tmp_a, *coeff);
346        mpq_mul(*tmp_b, *tmp_a, *tmp_b);
347      }
348
349      mpq_add(*mp_constant, *mp_constant, *tmp_b);     
350    } else {
351      fprintf(stdout, "ERROR: WRONG LINEAR ARITHMETIC ATOM TYPE\n");
352      str_a = smt_convert_fml_to_string(lfml);
353      str_b = smt_convert_fml_to_string(rfml);
354      fprintf(stdout, "fml = %s * %s\n", str_a, str_b);
355      fprintf(stdout, "unknown\n");
356      exit(0);     
357    }
358  }
359 
360  return;
361}
362
363int
364smt_mp_dl_theory_solve(smtManager_t * sm)
365{
366  smt_mp_bellman_ford_main(sm);
367
368  if (!(sm->flag & SAT_MASK)) { 
369
370    return 2; /* theory unsat */
371
372  }
373
374  if (sm->litArr->size == sm->avarArr->size) {
375
376    smt_mp_check_solution(sm);
377   
378  }
379   
380  return 0; /* sat without theory prop */
381}
382
383void
384smt_mp_bellman_ford_main(smtManager_t * sm)
385{
386  sm->stats->num_bf_call++;
387
388  smt_mp_generate_constraint_graph(sm);
389
390  smt_mp_bellman_ford_algorithm(sm);
391}
392
393void
394smt_mp_generate_constraint_graph(smtManager_t * sm)
395{
396  smtMp_t * mp = sm->mp;
397  mpq_t * pool = mp->pool;
398  mpq_t * weight = 0, * epsilon = 0;
399  smtGraph_t * G;
400  smtVertex_t * src, * dest;
401  smtEdge_t * e = 0;
402  smtAvar_t * avar;
403  smtNvar_t * lnvar, *rnvar;
404  satArray_t *cur_edges; 
405  int cur_index, lit, id;
406  int vindex, eindex;
407  int i;
408
409  /* mp vars */
410  weight = &pool[0];
411  epsilon = &pool[1];
412 
413  if (sm->cG == 0) smt_init_constraint_graph(sm);   
414
415  G = sm->cG;
416  eindex = G->esize;
417  cur_edges = G->cur_edges;
418  cur_edges->size = 0;
419
420  /* generate edges */
421  cur_index = sm->cur_index;
422
423  for(i = cur_index; i < sm->litArr->size; i++) {
424    lit = sm->litArr->space[i];
425    id = (lit>0) ? lit : -lit;
426    avar = (smtAvar_t *) sm->id2var->space[id];
427    /* filter theory propagated avar */
428    if (avar->flag & TPROP_MASK) continue;
429    assert(avar->type != EQ_c);
430     
431    if (sm->avalues[id] == 1) {
432      lnvar = (smtNvar_t *) avar->nvars->space[0];
433      rnvar = (smtNvar_t *) avar->nvars->space[1];
434      mpq_set(*weight, mp->constants[avar->id]);
435    } else {
436      assert(sm->avalues[id] == 0);
437      lnvar = (smtNvar_t *) avar->nvars->space[1];
438      rnvar = (smtNvar_t *) avar->nvars->space[0];
439      mpq_neg(*weight, mp->constants[e->avar->id]);
440      mpq_set_d(*epsilon, sm->epsilon);
441      mpq_sub(*weight, *weight, *epsilon);
442    }
443
444    vindex = lnvar->id - 1;
445    dest = &(G->vHead[vindex]);
446    vindex = rnvar->id - 1;
447    src =  &(G->vHead[vindex]);
448
449    e = smt_find_edge(src, dest);
450
451    if (e) {
452      if ( mpq_cmp(*weight, mp->weights[avar->id]) < 0 ) {
453        e->avar->flag |= IMPLIED_MASK;
454        e->implied = sat_array_insert(e->implied, (long) e->avar);
455        e->avar = avar;
456        mpq_set(mp->weights[avar->id], *weight);
457        cur_edges = sat_array_insert(cur_edges, (long) e);
458
459      } else {
460        avar->flag |= IMPLIED_MASK;
461        e->implied = sat_array_insert(e->implied, (long) avar);
462      }
463    } else {
464      e = smt_add_edge(G, src, dest, avar, eindex++);
465      mpq_set(mp->weights[avar->id], *weight);
466      cur_edges = sat_array_insert(cur_edges, (long) e);
467      G->esize++;
468    }
469  }
470
471  return;
472}
473
474void
475smt_mp_bellman_ford_algorithm(smtManager_t * sm)     
476{
477  smtMp_t * mp = sm->mp;
478  mpq_t * pool = mp->pool;
479  mpq_t * new_dist = 0, * cur_dist = 0, * weight = 0;
480  smtGraph_t * G = sm->cG;
481  smtVertex_t * v = 0, * src = 0, * dest = 0, *parent = 0;
482  smtEdge_t * e;
483  smtQueue_t * Q;
484  int i, qsize, cycleFound = 0;
485
486  /* mp vars */
487  new_dist = &pool[0];
488  cur_dist = &pool[1];
489  weight = &pool[2];
490 
491  smt_mp_init_bellman_ford_algorithm(sm); 
492
493  Q = G->queue;
494
495  qsize = Q->size;
496
497  while( (v = (smtVertex_t *) smt_dequeue(Q)) ) {
498
499    G->flags[v->index] &= RESET_FRINGE_MASK;
500
501    if (G->flags[v->index] & VISITED_MASK)
502      continue;
503
504    for(i = 0; i < v->outs->size; i++) {
505      e = (smtEdge_t *) v->outs->space[i];
506      src = e->src;
507      dest = e->dest;
508
509      mpq_set(*weight, mp->weights[e->avar->id]);
510      mpq_add(*new_dist, mp->values[src->index], *weight);
511      mpq_set(*cur_dist, mp->values[dest->index]);
512
513      if ( mpq_cmp(*new_dist, *cur_dist) < 0 ) {
514        /* check if src is in the subtree of dest.
515           if this is the case, the negative cycle is detected */
516        parent = G->paths[v->index];
517        while(1) {
518          if (parent == 0)
519            break;
520          else if(parent == dest) {
521            G->paths[dest->index] = v;
522            G->flags[dest->index] &= RESET_VISITED_MASK;
523            cycleFound = 1;
524            break;
525          }
526          parent = G->paths[parent->index];
527        }
528        if (cycleFound) break;
529
530        smt_delete_subtree(G, dest);
531
532        /* relaxation */
533        mpq_set(mp->values[dest->index], *new_dist);
534       
535        G->paths[dest->index] = src;
536        G->flags[dest->index] &= RESET_VISITED_MASK;
537       
538        if (!(G->flags[dest->index] & FRINGE_MASK)) {
539          G->flags[dest->index] |= FRINGE_MASK; 
540          smt_enqueue(Q, (long) dest);
541        }
542      }
543    }
544
545    if(cycleFound)      break;
546  }
547
548  if(cycleFound) {
549
550    sm->stats->num_bf_conf++;
551
552    sm->flag &= RESET_SAT_MASK;
553
554    smt_collect_edges_in_neg_cycle(G, dest);
555
556    smt_get_lemma_from_neg_cycle(sm, G->neg_edges);
557
558    smt_mp_retrieve_previous_distance(sm);
559
560  } else if (qsize) {
561   
562    smt_update_value_with_current_distance(sm);
563
564  }
565
566  return;
567}
568
569void
570smt_mp_init_bellman_ford_algorithm(smtManager_t * sm)
571{
572  smtMp_t * mp = sm->mp;
573  mpq_t * pool = mp->pool;
574  mpq_t * new_dist = 0, * cur_dist = 0, * weight = 0;
575  smtGraph_t * G = sm->cG;
576  smtVertex_t * src, * dest;
577  smtEdge_t * e;
578  smtQueue_t * Q;
579  int i;
580
581  /* mp vars */
582  new_dist = &pool[0];
583  cur_dist = &pool[1];
584  weight = &pool[2];
585 
586  memset(G->flags, 0, sizeof(int) * G->vsize);
587  memset(G->paths, 0, sizeof(smtVertex_t *) * G->vsize);
588
589  Q = G->queue;
590  smt_init_queue(Q);
591
592  for(i = 0; i < G->cur_edges->size; i++) {
593    e = (smtEdge_t *) G->cur_edges->space[i];
594    src = e->src;
595    dest = e->dest;
596
597    mpq_set(*weight, mp->weights[e->avar->id]);
598    mpq_add(*new_dist, mp->values[src->index], *weight);
599    mpq_set(*cur_dist, mp->values[dest->index]);
600
601    if ( mpq_cmp(*new_dist, *cur_dist) < 0 ) {
602      if(!(G->flags[src->index] & FRINGE_MASK)) {
603        G->flags[src->index] |= FRINGE_MASK; /* not in queue */
604        smt_enqueue(Q, (long) src);
605      }
606    }
607  }
608
609  return;
610}
611
612void
613smt_mp_retrieve_previous_distance(smtManager_t * sm)
614{
615  smtMp_t * mp = sm->mp;
616  int i, size;
617
618  size = sm->nvarArr->size + 1;
619
620  for(i = 0; i < size; i++) {
621    mpq_set(mp->pvalues[i], mp->values[i]);
622  }
623
624  return;
625}
626
627void
628smt_mp_print_atom_constants(smtManager_t * sm)
629{
630  smtMp_t * mp = sm->mp;
631  smtAvar_t * avar;
632  int i, id;
633
634  for(i = 0; i < sm->avarArr->size; i++) {
635    avar = (smtAvar_t *) sm->avarArr->space[i];
636    id = avar->id;
637    fprintf(stdout, "%s\n", avar->name);
638    gmp_printf("constant:%#Qd\n", mp->constants[id]);
639  }
640}
641
642void
643smt_mp_print_value(mpq_t * value)
644{
645  gmp_printf("value:%#Qd\n", *value); 
646}
647#endif
Note: See TracBrowser for help on using the repository browser.