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
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
|
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html><head><meta http-equiv="Content-Type" content="text/html;charset=UTF-8">
<title>polylib: alpha.c Source File</title>
<link href="doxygen.css" rel="stylesheet" type="text/css">
<link href="tabs.css" rel="stylesheet" type="text/css">
</head><body>
<!-- Generated by Doxygen 1.5.6 -->
<div class="navigation" id="top">
<div class="tabs">
<ul>
<li><a href="main.html"><span>Main Page</span></a></li>
<li><a href="annotated.html"><span>Classes</span></a></li>
<li class="current"><a href="files.html"><span>Files</span></a></li>
</ul>
</div>
<h1>alpha.c</h1><a href="alpha_8c.html">Go to the documentation of this file.</a><div class="fragment"><pre class="fragment"><a name="l00001"></a>00001 <span class="comment">/* polytest.c */</span>
<a name="l00002"></a>00002 <span class="preprocessor">#include <stdio.h></span>
<a name="l00003"></a>00003 <span class="preprocessor">#include <<a class="code" href="polylib_8h.html">polylib/polylib.h</a>></span>
<a name="l00004"></a>00004
<a name="l00005"></a>00005 <span class="comment">/* alpha.c</span>
<a name="l00006"></a>00006 <span class="comment"> COPYRIGHT</span>
<a name="l00007"></a>00007 <span class="comment"> Both this software and its documentation are</span>
<a name="l00008"></a>00008 <span class="comment"></span>
<a name="l00009"></a>00009 <span class="comment"> Copyright 1993 by IRISA /Universite de Rennes I - France,</span>
<a name="l00010"></a>00010 <span class="comment"> Copyright 1995,1996 by BYU, Provo, Utah</span>
<a name="l00011"></a>00011 <span class="comment"> all rights reserved.</span>
<a name="l00012"></a>00012 <span class="comment"></span>
<a name="l00013"></a>00013 <span class="comment"> Permission is granted to copy, use, and distribute</span>
<a name="l00014"></a>00014 <span class="comment"> for any commercial or noncommercial purpose under the terms</span>
<a name="l00015"></a>00015 <span class="comment"> of the GNU General Public license, version 2, June 1991</span>
<a name="l00016"></a>00016 <span class="comment"> (see file : LICENSING).</span>
<a name="l00017"></a>00017 <span class="comment">*/</span>
<a name="l00018"></a>00018
<a name="l00019"></a>00019 <span class="preprocessor">#include <stdio.h></span>
<a name="l00020"></a>00020 <span class="preprocessor">#include <stdlib.h></span>
<a name="l00021"></a>00021 <span class="preprocessor">#include <string.h></span>
<a name="l00022"></a>00022
<a name="l00023"></a>00023 <span class="comment">/*---------------------------------------------------------------------*/</span>
<a name="l00024"></a>00024 <span class="comment">/* int exist_points(pos,P,context) */</span>
<a name="l00025"></a>00025 <span class="comment">/* pos : index position of current loop index (0..hdim-1) */</span>
<a name="l00026"></a>00026 <span class="comment">/* P: loop domain */</span>
<a name="l00027"></a>00027 <span class="comment">/* context : context values for fixed indices */</span>
<a name="l00028"></a>00028 <span class="comment">/* recursive procedure, recurs for each imbriquation */</span>
<a name="l00029"></a>00029 <span class="comment">/* returns 1 if there exists any integer points in this polyhedron */</span>
<a name="l00030"></a>00030 <span class="comment">/* returns 0 if there are none */</span>
<a name="l00031"></a>00031 <span class="comment">/*---------------------------------------------------------------------*/</span>
<a name="l00032"></a><a class="code" href="alpha_8c.html#3a7bf711585ac2e71c8cfab9e3af672b">00032</a> <span class="keyword">static</span> <span class="keywordtype">int</span> <a class="code" href="alpha_8c.html#3a7bf711585ac2e71c8cfab9e3af672b">exist_points</a>(<span class="keywordtype">int</span> pos,<a class="code" href="structpolyhedron.html">Polyhedron</a> *Pol,Value *context) {
<a name="l00033"></a>00033
<a name="l00034"></a>00034 Value LB, UB, k,tmp;
<a name="l00035"></a>00035
<a name="l00036"></a>00036 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(LB); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(UB);
<a name="l00037"></a>00037 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(k); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(tmp);
<a name="l00038"></a>00038 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(LB,0);
<a name="l00039"></a>00039 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(UB,0);
<a name="l00040"></a>00040
<a name="l00041"></a>00041 <span class="comment">/* Problem if UB or LB is INFINITY */</span>
<a name="l00042"></a>00042 <span class="keywordflow">if</span> (<a class="code" href="polyhedron_8c.html#a8c3f429dd74c5ec2087f97d702cdc66">lower_upper_bounds</a>(pos,Pol,context,&LB,&UB) !=0) {
<a name="l00043"></a>00043 <a class="code" href="errormsg_8c.html#1ac398d9106568a1991a001172e0e00b">errormsg1</a>(<span class="stringliteral">"exist_points"</span>, <span class="stringliteral">"infdom"</span>, <span class="stringliteral">"infinite domain"</span>);
<a name="l00044"></a>00044 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(LB);
<a name="l00045"></a>00045 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(UB);
<a name="l00046"></a>00046 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(k);
<a name="l00047"></a>00047 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(tmp);
<a name="l00048"></a>00048 <span class="keywordflow">return</span> -1;
<a name="l00049"></a>00049 }
<a name="l00050"></a>00050 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(context[pos],0);
<a name="l00051"></a>00051 <span class="keywordflow">if</span>(<a class="code" href="source_2arith_2arithmetique_8h.html#d87a1fcce82a9885ffd8518f4fba2c0e">value_lt</a>(UB,LB)) {
<a name="l00052"></a>00052 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(LB);
<a name="l00053"></a>00053 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(UB);
<a name="l00054"></a>00054 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(k);
<a name="l00055"></a>00055 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(tmp);
<a name="l00056"></a>00056 <span class="keywordflow">return</span> 0;
<a name="l00057"></a>00057 }
<a name="l00058"></a>00058 <span class="keywordflow">if</span> (!Pol-><a class="code" href="structpolyhedron.html#ab97126cb9e56431d40eaf4e68b3953d">next</a>) {
<a name="l00059"></a>00059 <a class="code" href="source_2arith_2arithmetique_8h.html#e92a58eee3b6f5c6a99e6837e68407e1">value_subtract</a>(tmp,UB,LB);
<a name="l00060"></a>00060 <a class="code" href="source_2arith_2arithmetique_8h.html#88693f35dd41deddc6ac0700073fc8db">value_increment</a>(tmp,tmp);
<a name="l00061"></a>00061 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(UB);
<a name="l00062"></a>00062 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(LB);
<a name="l00063"></a>00063 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(k);
<a name="l00064"></a>00064 <span class="keywordflow">return</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#01904c8c75c10d407a24c55f16ec27c7">value_pos_p</a>(tmp));
<a name="l00065"></a>00065 }
<a name="l00066"></a>00066
<a name="l00067"></a>00067 <span class="keywordflow">for</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(k,LB);<a class="code" href="source_2arith_2arithmetique_8h.html#47975ace017981602e1064f98f43f8a7">value_le</a>(k,UB);<a class="code" href="source_2arith_2arithmetique_8h.html#88693f35dd41deddc6ac0700073fc8db">value_increment</a>(k,k)) {
<a name="l00068"></a>00068
<a name="l00069"></a>00069 <span class="comment">/* insert k in context */</span>
<a name="l00070"></a>00070 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(context[pos],k);
<a name="l00071"></a>00071 <span class="keywordflow">if</span> (<a class="code" href="alpha_8c.html#3a7bf711585ac2e71c8cfab9e3af672b">exist_points</a>(pos+1,Pol-><a class="code" href="structpolyhedron.html#ab97126cb9e56431d40eaf4e68b3953d">next</a>,context) > 0 ) {
<a name="l00072"></a>00072 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(LB); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(UB);
<a name="l00073"></a>00073 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(k); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(tmp);
<a name="l00074"></a>00074 <span class="keywordflow">return</span> 1;
<a name="l00075"></a>00075 }
<a name="l00076"></a>00076 }
<a name="l00077"></a>00077 <span class="comment">/* Reset context */</span>
<a name="l00078"></a>00078 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(context[pos],0);
<a name="l00079"></a>00079 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(UB); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(LB);
<a name="l00080"></a>00080 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(k); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(tmp);
<a name="l00081"></a>00081 <span class="keywordflow">return</span> 0;
<a name="l00082"></a>00082 }
<a name="l00083"></a>00083
<a name="l00084"></a>00084 <span class="comment">/*--------------------------------------------------------------*/</span>
<a name="l00085"></a>00085 <span class="comment">/* Test to see if there are any integral points in a polyhedron */</span>
<a name="l00086"></a>00086 <span class="comment">/*--------------------------------------------------------------*/</span>
<a name="l00087"></a><a class="code" href="alpha_8h.html#5be7eb3e254a2e2e61c50d736b4900de">00087</a> <span class="keywordtype">int</span> <a class="code" href="alpha_8c.html#5be7eb3e254a2e2e61c50d736b4900de">Polyhedron_Not_Empty</a>(<a class="code" href="structpolyhedron.html">Polyhedron</a> *P,<a class="code" href="structpolyhedron.html">Polyhedron</a> *C,<span class="keywordtype">int</span> <a class="code" href="verif__ehrhart_8c.html#89fd83aa168651629c012d8655635588">MAXRAYS</a>) {
<a name="l00088"></a>00088
<a name="l00089"></a>00089 <span class="keywordtype">int</span> res,i;
<a name="l00090"></a>00090 Value *context;
<a name="l00091"></a>00091 <a class="code" href="structpolyhedron.html">Polyhedron</a> *L;
<a name="l00092"></a>00092
<a name="l00093"></a>00093 <a class="code" href="polyhedron_8h.html#729f721c450a034c1219447745c6b123">POL_ENSURE_FACETS</a>(P);
<a name="l00094"></a>00094 <a class="code" href="polyhedron_8h.html#317e5eb888f5c14e1229742f31169378">POL_ENSURE_VERTICES</a>(P);
<a name="l00095"></a>00095 <a class="code" href="polyhedron_8h.html#729f721c450a034c1219447745c6b123">POL_ENSURE_FACETS</a>(C);
<a name="l00096"></a>00096 <a class="code" href="polyhedron_8h.html#317e5eb888f5c14e1229742f31169378">POL_ENSURE_VERTICES</a>(C);
<a name="l00097"></a>00097
<a name="l00098"></a>00098 <span class="comment">/* Create a context vector size dim+2 and set it to all zeros */</span>
<a name="l00099"></a>00099 context = (Value *) malloc((P-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>+2)*<span class="keyword">sizeof</span>(Value));
<a name="l00100"></a>00100
<a name="l00101"></a>00101 <span class="comment">/* Initialize array 'context' */</span>
<a name="l00102"></a>00102 <span class="keywordflow">for</span> (i=0;i<(P-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>+2);i++)
<a name="l00103"></a>00103 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(context[i]);
<a name="l00104"></a>00104
<a name="l00105"></a>00105 <a class="code" href="vector_8c.html#27ccb3ea01f3a496bb799fb99e9e3075">Vector_Set</a>(context,0,(P-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>+2));
<a name="l00106"></a>00106
<a name="l00107"></a>00107 <span class="comment">/* Set context[P->Dimension+1] = 1 (the constant) */</span>
<a name="l00108"></a>00108 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(context[P-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>+1],1);
<a name="l00109"></a>00109
<a name="l00110"></a>00110 L = <a class="code" href="polyhedron_8c.html#7e6c09758e3d3063be1386c7026b38ba">Polyhedron_Scan</a>(P,C,MAXRAYS);
<a name="l00111"></a>00111 res = <a class="code" href="alpha_8c.html#3a7bf711585ac2e71c8cfab9e3af672b">exist_points</a>(1,L,context);
<a name="l00112"></a>00112 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(L);
<a name="l00113"></a>00113
<a name="l00114"></a>00114 <span class="comment">/* Clear array 'context' */</span>
<a name="l00115"></a>00115 <span class="keywordflow">for</span> (i=0;i<(P-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>+2);i++)
<a name="l00116"></a>00116 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(context[i]);
<a name="l00117"></a>00117 free(context);
<a name="l00118"></a>00118 <span class="keywordflow">return</span> res;
<a name="l00119"></a>00119 }
<a name="l00120"></a>00120
<a name="l00121"></a>00121 <span class="comment">/* PolyhedronLTQ ( P1, P2 ) */</span>
<a name="l00122"></a>00122 <span class="comment">/* P1 and P2 must be simple polyhedra */</span>
<a name="l00123"></a>00123 <span class="comment">/* result = 0 --> not comparable */</span>
<a name="l00124"></a>00124 <span class="comment">/* result = -1 --> P1 < P2 */</span>
<a name="l00125"></a>00125 <span class="comment">/* result = 1 --> P1 > P2 */</span>
<a name="l00126"></a>00126 <span class="comment">/* INDEX = 1 .... Dimension */</span>
<a name="l00127"></a><a class="code" href="alpha_8h.html#f6fa8b8fe252ab66fdd61f14595ff979">00127</a> <span class="keywordtype">int</span> <a class="code" href="alpha_8c.html#0f28b3a3172da3ef981f4a57c7ead628">PolyhedronLTQ</a> (<a class="code" href="structpolyhedron.html">Polyhedron</a> *Pol1,<a class="code" href="structpolyhedron.html">Polyhedron</a> *Pol2,<span class="keywordtype">int</span> INDEX, <span class="keywordtype">int</span> PDIM, <span class="keywordtype">int</span> NbMaxConstrs) {
<a name="l00128"></a>00128
<a name="l00129"></a>00129 <span class="keywordtype">int</span> res, dim, i, j, k;
<a name="l00130"></a>00130 <a class="code" href="structpolyhedron.html">Polyhedron</a> *Q1, *Q2, *Q3, *Q4, *Q;
<a name="l00131"></a>00131 <a class="code" href="structmatrix.html">Matrix</a> *Mat;
<a name="l00132"></a>00132
<a name="l00133"></a>00133 <span class="keywordflow">if</span> (Pol1-><a class="code" href="structpolyhedron.html#ab97126cb9e56431d40eaf4e68b3953d">next</a> || Pol2-><a class="code" href="structpolyhedron.html#ab97126cb9e56431d40eaf4e68b3953d">next</a>) {
<a name="l00134"></a>00134 <a class="code" href="errormsg_8c.html#1ac398d9106568a1991a001172e0e00b">errormsg1</a>(<span class="stringliteral">"PolyhedronLTQ"</span>, <span class="stringliteral">"compoly"</span>, <span class="stringliteral">"Can only compare polyhedra"</span>);
<a name="l00135"></a>00135 <span class="keywordflow">return</span> 0;
<a name="l00136"></a>00136 }
<a name="l00137"></a>00137 <span class="keywordflow">if</span> (Pol1-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a> != Pol2-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>) {
<a name="l00138"></a>00138 <a class="code" href="errormsg_8c.html#1ac398d9106568a1991a001172e0e00b">errormsg1</a>(<span class="stringliteral">"PolyhedronLTQ"</span>,<span class="stringliteral">"diffdim"</span>,<span class="stringliteral">"Polyhedra are not same dimension"</span>);
<a name="l00139"></a>00139 <span class="keywordflow">return</span> 0;
<a name="l00140"></a>00140 }
<a name="l00141"></a>00141 dim = Pol1-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>+2;
<a name="l00142"></a>00142
<a name="l00143"></a>00143 <a class="code" href="polyhedron_8h.html#729f721c450a034c1219447745c6b123">POL_ENSURE_FACETS</a>(Pol1);
<a name="l00144"></a>00144 <a class="code" href="polyhedron_8h.html#317e5eb888f5c14e1229742f31169378">POL_ENSURE_VERTICES</a>(Pol1);
<a name="l00145"></a>00145 <a class="code" href="polyhedron_8h.html#729f721c450a034c1219447745c6b123">POL_ENSURE_FACETS</a>(Pol2);
<a name="l00146"></a>00146 <a class="code" href="polyhedron_8h.html#317e5eb888f5c14e1229742f31169378">POL_ENSURE_VERTICES</a>(Pol2);
<a name="l00147"></a>00147
<a name="l00148"></a>00148 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00149"></a>00149 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"P1\n"</span>);
<a name="l00150"></a>00150 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Pol1);
<a name="l00151"></a>00151 fprintf(stdout, <span class="stringliteral">"P2\n"</span>);
<a name="l00152"></a>00152 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Pol2);
<a name="l00153"></a>00153 <span class="preprocessor">#endif</span>
<a name="l00154"></a>00154 <span class="preprocessor"></span>
<a name="l00155"></a>00155 <span class="comment">/* Create the Line to add */</span>
<a name="l00156"></a>00156 k = Pol1-><a class="code" href="structpolyhedron.html#2a02cea8b7ba3dde415041b8b2373bc8">Dimension</a>-INDEX+1-PDIM;
<a name="l00157"></a>00157 Mat = <a class="code" href="matrix_8c.html#c0b29e1d99a2823ad00b5f2157879d80">Matrix_Alloc</a>(k,dim);
<a name="l00158"></a>00158 <a class="code" href="vector_8c.html#27ccb3ea01f3a496bb799fb99e9e3075">Vector_Set</a>(Mat-><a class="code" href="structmatrix.html#9c294a3fd4d273963e27a0c8cd819bd5">p_Init</a>,0,dim*k);
<a name="l00159"></a>00159 <span class="keywordflow">for</span>(j=0,i=INDEX;j<k;i++,j++)
<a name="l00160"></a>00160 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[j][i],1);
<a name="l00161"></a>00161
<a name="l00162"></a>00162 Q1 = <a class="code" href="polyhedron_8c.html#964e678fe2d07dd6162a7f1c429d38ba">AddRays</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0],k,Pol1,NbMaxConstrs);
<a name="l00163"></a>00163 Q2 = <a class="code" href="polyhedron_8c.html#964e678fe2d07dd6162a7f1c429d38ba">AddRays</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0],k,Pol2,NbMaxConstrs);
<a name="l00164"></a>00164
<a name="l00165"></a>00165 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00166"></a>00166 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"Q1\n"</span>);
<a name="l00167"></a>00167 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Q1);
<a name="l00168"></a>00168 fprintf(stdout, <span class="stringliteral">"Q2\n"</span>);
<a name="l00169"></a>00169 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Q2);
<a name="l00170"></a>00170 <span class="preprocessor">#endif</span>
<a name="l00171"></a>00171 <span class="preprocessor"></span>
<a name="l00172"></a>00172 <a class="code" href="matrix_8c.html#fcb312b7c12a6997cd66964ecc34e1a6">Matrix_Free</a>(Mat);
<a name="l00173"></a>00173 Q = <a class="code" href="polyhedron_8c.html#c5a1a2751f0b833183560af25b18033b">DomainIntersection</a>(Q1,Q2,NbMaxConstrs);
<a name="l00174"></a>00174
<a name="l00175"></a>00175 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00176"></a>00176 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"Q\n"</span>);
<a name="l00177"></a>00177 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Q);
<a name="l00178"></a>00178 <span class="preprocessor">#endif</span>
<a name="l00179"></a>00179 <span class="preprocessor"></span>
<a name="l00180"></a>00180 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q1);
<a name="l00181"></a>00181 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q2);
<a name="l00182"></a>00182
<a name="l00183"></a>00183 <span class="keywordflow">if</span> (<a class="code" href="types_8h.html#156641ae25106f9a6ac5c5c2a624c168">emptyQ</a>(Q)) res = 0; <span class="comment">/* not comparable */</span>
<a name="l00184"></a>00184 <span class="keywordflow">else</span> {
<a name="l00185"></a>00185 Q1 = <a class="code" href="polyhedron_8c.html#c5a1a2751f0b833183560af25b18033b">DomainIntersection</a>(Pol1,Q,NbMaxConstrs);
<a name="l00186"></a>00186 Q2 = <a class="code" href="polyhedron_8c.html#c5a1a2751f0b833183560af25b18033b">DomainIntersection</a>(Pol2,Q,NbMaxConstrs);
<a name="l00187"></a>00187
<a name="l00188"></a>00188 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00189"></a>00189 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"Q1\n"</span>);
<a name="l00190"></a>00190 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Q1);
<a name="l00191"></a>00191 fprintf(stdout, <span class="stringliteral">"Q2\n"</span>);
<a name="l00192"></a>00192 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Q2);
<a name="l00193"></a>00193 <span class="preprocessor">#endif</span>
<a name="l00194"></a>00194 <span class="preprocessor"></span>
<a name="l00195"></a>00195 k = Q1-><a class="code" href="structpolyhedron.html#2db05293e731089b9198d45cf027bd1c">NbConstraints</a> + Q2-><a class="code" href="structpolyhedron.html#2db05293e731089b9198d45cf027bd1c">NbConstraints</a>;
<a name="l00196"></a>00196 Mat = <a class="code" href="matrix_8c.html#c0b29e1d99a2823ad00b5f2157879d80">Matrix_Alloc</a>(k, dim);
<a name="l00197"></a>00197 <a class="code" href="vector_8c.html#27ccb3ea01f3a496bb799fb99e9e3075">Vector_Set</a>(Mat-><a class="code" href="structmatrix.html#9c294a3fd4d273963e27a0c8cd819bd5">p_Init</a>,0,k*dim);
<a name="l00198"></a>00198
<a name="l00199"></a>00199 <span class="comment">/* First compute surrounding polyhedron */</span>
<a name="l00200"></a>00200 j=0;
<a name="l00201"></a>00201 <span class="keywordflow">for</span> (i=0; i<Q1-><a class="code" href="structpolyhedron.html#2db05293e731089b9198d45cf027bd1c">NbConstraints</a>; i++) {
<a name="l00202"></a>00202 <span class="keywordflow">if</span> ((<a class="code" href="source_2arith_2arithmetique_8h.html#61db0ba9ba925615ec70eb92e40b9ef4">value_one_p</a>(Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][0])) && (<a class="code" href="source_2arith_2arithmetique_8h.html#01904c8c75c10d407a24c55f16ec27c7">value_pos_p</a>(Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][INDEX]))) {
<a name="l00203"></a>00203
<a name="l00204"></a>00204 <span class="comment">/* keep Q1's lower bounds */</span>
<a name="l00205"></a>00205 <span class="keywordflow">for</span> (k=0; k<dim; k++)
<a name="l00206"></a>00206 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[j][k],Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][k]);
<a name="l00207"></a>00207 j++;
<a name="l00208"></a>00208 }
<a name="l00209"></a>00209 }
<a name="l00210"></a>00210 <span class="keywordflow">for</span> (i=0; i<Q2-><a class="code" href="structpolyhedron.html#2db05293e731089b9198d45cf027bd1c">NbConstraints</a>; i++) {
<a name="l00211"></a>00211 <span class="keywordflow">if</span> ((<a class="code" href="source_2arith_2arithmetique_8h.html#61db0ba9ba925615ec70eb92e40b9ef4">value_one_p</a>(Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][0])) && (<a class="code" href="source_2arith_2arithmetique_8h.html#f031da52c8a160fd1ee6e4c793f8c78a">value_neg_p</a>(Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][INDEX]))) {
<a name="l00212"></a>00212
<a name="l00213"></a>00213 <span class="comment">/* and keep Q2's upper bounds */</span>
<a name="l00214"></a>00214 <span class="keywordflow">for</span> (k=0; k<dim; k++)
<a name="l00215"></a>00215 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[j][k],Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][k]);
<a name="l00216"></a>00216 j++;
<a name="l00217"></a>00217 }
<a name="l00218"></a>00218 }
<a name="l00219"></a>00219 Q4 = <a class="code" href="polyhedron_8c.html#9690baafb863abbd2e8d5c4da2995416">AddConstraints</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0], j, Q, NbMaxConstrs);
<a name="l00220"></a>00220 <a class="code" href="matrix_8c.html#fcb312b7c12a6997cd66964ecc34e1a6">Matrix_Free</a>(Mat);
<a name="l00221"></a>00221
<a name="l00222"></a>00222 <span class="preprocessor">#ifdef debug</span>
<a name="l00223"></a>00223 <span class="preprocessor"></span> fprintf(stderr, <span class="stringliteral">"Q4 surrounding polyhedron\n"</span>);
<a name="l00224"></a>00224 Polyhderon_Print(stderr,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>, Q4);
<a name="l00225"></a>00225 <span class="preprocessor">#endif</span>
<a name="l00226"></a>00226 <span class="preprocessor"></span>
<a name="l00227"></a>00227 <span class="comment">/* if surrounding polyhedron is empty, D1>D2 */</span>
<a name="l00228"></a>00228 <span class="keywordflow">if</span> (<a class="code" href="types_8h.html#156641ae25106f9a6ac5c5c2a624c168">emptyQ</a>(Q4)) {
<a name="l00229"></a>00229 res = 1;
<a name="l00230"></a>00230
<a name="l00231"></a>00231 <span class="preprocessor">#ifdef debug</span>
<a name="l00232"></a>00232 <span class="preprocessor"></span> fprintf(stderr, <span class="stringliteral">"Surrounding polyhedron is empty\n"</span>);
<a name="l00233"></a>00233 <span class="preprocessor">#endif</span>
<a name="l00234"></a>00234 <span class="preprocessor"></span> <span class="keywordflow">goto</span> LTQdone2;
<a name="l00235"></a>00235 }
<a name="l00236"></a>00236
<a name="l00237"></a>00237 <span class="comment">/* Test if Q1 < Q2 */</span>
<a name="l00238"></a>00238 <span class="comment">/* Build a constraint array for >= Q1 and <= Q2 */</span>
<a name="l00239"></a>00239 Mat = <a class="code" href="matrix_8c.html#c0b29e1d99a2823ad00b5f2157879d80">Matrix_Alloc</a>(2,dim);
<a name="l00240"></a>00240 <a class="code" href="vector_8c.html#27ccb3ea01f3a496bb799fb99e9e3075">Vector_Set</a>(Mat-><a class="code" href="structmatrix.html#9c294a3fd4d273963e27a0c8cd819bd5">p_Init</a>,0,2*dim);
<a name="l00241"></a>00241
<a name="l00242"></a>00242 <span class="comment">/* Choose a contraint from Q1 */</span>
<a name="l00243"></a>00243 <span class="keywordflow">for</span> (i=0; i<Q1-><a class="code" href="structpolyhedron.html#2db05293e731089b9198d45cf027bd1c">NbConstraints</a>; i++) {
<a name="l00244"></a>00244 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#827532f2140ae2aa96e46baebae09723">value_zero_p</a>(Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][0])) {
<a name="l00245"></a>00245
<a name="l00246"></a>00246 <span class="comment">/* Equality */</span>
<a name="l00247"></a>00247 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#827532f2140ae2aa96e46baebae09723">value_zero_p</a>(Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][INDEX])) {
<a name="l00248"></a>00248
<a name="l00249"></a>00249 <span class="comment">/* Ignore side constraint (they are in Q) */</span>
<a name="l00250"></a>00250 <span class="keywordflow">continue</span>;
<a name="l00251"></a>00251 }
<a name="l00252"></a>00252 <span class="keywordflow">else</span> <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#f031da52c8a160fd1ee6e4c793f8c78a">value_neg_p</a>(Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][INDEX])) {
<a name="l00253"></a>00253
<a name="l00254"></a>00254 <span class="comment">/* copy -constraint to Mat */</span>
<a name="l00255"></a>00255 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0][0],1);
<a name="l00256"></a>00256 <span class="keywordflow">for</span> (k=1; k<dim; k++)
<a name="l00257"></a>00257 <a class="code" href="source_2arith_2arithmetique_8h.html#0daf9a8ecdc14e7a274832b0d2add830">value_oppose</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0][k],Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][k]);
<a name="l00258"></a>00258 }
<a name="l00259"></a>00259 <span class="keywordflow">else</span> {
<a name="l00260"></a>00260
<a name="l00261"></a>00261 <span class="comment">/* Copy constraint to Mat */</span>
<a name="l00262"></a>00262
<a name="l00263"></a>00263 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0][0],1);
<a name="l00264"></a>00264 <span class="keywordflow">for</span> (k=1; k<dim; k++)
<a name="l00265"></a>00265 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0][k],Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][k]);
<a name="l00266"></a>00266 }
<a name="l00267"></a>00267 }
<a name="l00268"></a>00268 <span class="keywordflow">else</span> <span class="keywordflow">if</span>(<a class="code" href="source_2arith_2arithmetique_8h.html#f031da52c8a160fd1ee6e4c793f8c78a">value_neg_p</a>(Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][INDEX])) {
<a name="l00269"></a>00269
<a name="l00270"></a>00270 <span class="comment">/* Upper bound -- make a lower bound from it */</span>
<a name="l00271"></a>00271 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0][0],1);
<a name="l00272"></a>00272 <span class="keywordflow">for</span> (k=1; k<dim; k++)
<a name="l00273"></a>00273 <a class="code" href="source_2arith_2arithmetique_8h.html#0daf9a8ecdc14e7a274832b0d2add830">value_oppose</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0][k],Q1-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[i][k]);
<a name="l00274"></a>00274 }
<a name="l00275"></a>00275 <span class="keywordflow">else</span> {
<a name="l00276"></a>00276
<a name="l00277"></a>00277 <span class="comment">/* Lower or side bound -- ignore it */</span>
<a name="l00278"></a>00278 <span class="keywordflow">continue</span>;
<a name="l00279"></a>00279 }
<a name="l00280"></a>00280
<a name="l00281"></a>00281 <span class="comment">/* Choose a constraint from Q2 */</span>
<a name="l00282"></a>00282 <span class="keywordflow">for</span> (j=0; j<Q2-><a class="code" href="structpolyhedron.html#2db05293e731089b9198d45cf027bd1c">NbConstraints</a>; j++) {
<a name="l00283"></a>00283 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#827532f2140ae2aa96e46baebae09723">value_zero_p</a>(Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[j][0])) { <span class="comment">/* equality */</span>
<a name="l00284"></a>00284 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#827532f2140ae2aa96e46baebae09723">value_zero_p</a>(Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[j][INDEX])) {
<a name="l00285"></a>00285
<a name="l00286"></a>00286 <span class="comment">/* Ignore side constraint (they are in Q) */</span>
<a name="l00287"></a>00287 <span class="keywordflow">continue</span>;
<a name="l00288"></a>00288 }
<a name="l00289"></a>00289 <span class="keywordflow">else</span> <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#01904c8c75c10d407a24c55f16ec27c7">value_pos_p</a>(Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[j][INDEX])) {
<a name="l00290"></a>00290
<a name="l00291"></a>00291 <span class="comment">/* Copy -constraint to Mat */</span>
<a name="l00292"></a>00292 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[1][0],1);
<a name="l00293"></a>00293 <span class="keywordflow">for</span> (k=1; k<dim; k++)
<a name="l00294"></a>00294 <a class="code" href="source_2arith_2arithmetique_8h.html#0daf9a8ecdc14e7a274832b0d2add830">value_oppose</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[1][k],Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[j][k]);
<a name="l00295"></a>00295 }
<a name="l00296"></a>00296 <span class="keywordflow">else</span> {
<a name="l00297"></a>00297
<a name="l00298"></a>00298 <span class="comment">/* Copy constraint to Mat */</span>
<a name="l00299"></a>00299 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[1][0],1);
<a name="l00300"></a>00300 <span class="keywordflow">for</span> (k=1; k<dim; k++)
<a name="l00301"></a>00301 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[1][k],Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[j][k]);
<a name="l00302"></a>00302 };
<a name="l00303"></a>00303 }
<a name="l00304"></a>00304 <span class="keywordflow">else</span> <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#01904c8c75c10d407a24c55f16ec27c7">value_pos_p</a>(Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[j][INDEX])) {
<a name="l00305"></a>00305
<a name="l00306"></a>00306 <span class="comment">/* Lower bound -- make an upper bound from it */</span>
<a name="l00307"></a>00307 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[1][0],1);
<a name="l00308"></a>00308 <span class="keywordflow">for</span>(k=1;k<dim;k++)
<a name="l00309"></a>00309 <a class="code" href="source_2arith_2arithmetique_8h.html#0daf9a8ecdc14e7a274832b0d2add830">value_oppose</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[1][k],Q2-><a class="code" href="structpolyhedron.html#b2691b8eb5d2dc43a509ac2094d2bc73">Constraint</a>[j][k]);
<a name="l00310"></a>00310 }
<a name="l00311"></a>00311 <span class="keywordflow">else</span> {
<a name="l00312"></a>00312
<a name="l00313"></a>00313 <span class="comment">/* Upper or side bound -- ignore it */</span>
<a name="l00314"></a>00314 <span class="keywordflow">continue</span>;
<a name="l00315"></a>00315 };
<a name="l00316"></a>00316
<a name="l00317"></a>00317 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00318"></a>00318 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"i=%d j=%d M=\n"</span>, i+1, j+1);
<a name="l00319"></a>00319 <a class="code" href="matrix_8c.html#d4c3c350dba633f72bbcf3609b6b67d7">Matrix_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Mat);
<a name="l00320"></a>00320 <span class="preprocessor">#endif</span>
<a name="l00321"></a>00321 <span class="preprocessor"></span>
<a name="l00322"></a>00322 <span class="comment">/* Add Mat to Q and see if anything is made */</span>
<a name="l00323"></a>00323 Q3 = <a class="code" href="polyhedron_8c.html#9690baafb863abbd2e8d5c4da2995416">AddConstraints</a>(Mat-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[0],2,Q,NbMaxConstrs);
<a name="l00324"></a>00324
<a name="l00325"></a>00325 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00326"></a>00326 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"Q3\n"</span>);
<a name="l00327"></a>00327 <a class="code" href="polyhedron_8c.html#d87595bd01c5cdd0f3933824943db496">Polyhedron_Print</a>(stdout,<a class="code" href="types_8h.html#e6f16bcd4a42ba51cbb003e3d1e1cde6">P_VALUE_FMT</a>,Q3);
<a name="l00328"></a>00328 <span class="preprocessor">#endif</span>
<a name="l00329"></a>00329 <span class="preprocessor"></span>
<a name="l00330"></a>00330 <span class="keywordflow">if</span> (!<a class="code" href="types_8h.html#156641ae25106f9a6ac5c5c2a624c168">emptyQ</a>(Q3)) {
<a name="l00331"></a>00331 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q3);
<a name="l00332"></a>00332
<a name="l00333"></a>00333 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00334"></a>00334 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"not empty\n"</span>);
<a name="l00335"></a>00335 <span class="preprocessor">#endif</span>
<a name="l00336"></a>00336 <span class="preprocessor"></span> res = -1;
<a name="l00337"></a>00337 <span class="keywordflow">goto</span> LTQdone;
<a name="l00338"></a>00338 }
<a name="l00339"></a>00339 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00340"></a>00340 <span class="preprocessor"></span> fprintf(stdout,<span class="stringliteral">"empty\n"</span>);
<a name="l00341"></a>00341 <span class="preprocessor">#endif</span>
<a name="l00342"></a>00342 <span class="preprocessor"></span> <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q3);
<a name="l00343"></a>00343 } <span class="comment">/* end for j */</span>
<a name="l00344"></a>00344 } <span class="comment">/* end for i */</span>
<a name="l00345"></a>00345 res = 1;
<a name="l00346"></a>00346 LTQdone:
<a name="l00347"></a>00347 <a class="code" href="matrix_8c.html#fcb312b7c12a6997cd66964ecc34e1a6">Matrix_Free</a>(Mat);
<a name="l00348"></a>00348 LTQdone2:
<a name="l00349"></a>00349 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q4);
<a name="l00350"></a>00350 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q1);
<a name="l00351"></a>00351 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q2);
<a name="l00352"></a>00352 }
<a name="l00353"></a>00353 <a class="code" href="polyhedron_8c.html#e6d0a7daf8e801a777fc8e93d8cfe43a">Domain_Free</a>(Q);
<a name="l00354"></a>00354
<a name="l00355"></a>00355 <span class="preprocessor">#ifdef DEBUG</span>
<a name="l00356"></a>00356 <span class="preprocessor"></span> fprintf(stdout, <span class="stringliteral">"res = %d\n"</span>, res);
<a name="l00357"></a>00357 <span class="preprocessor">#endif</span>
<a name="l00358"></a>00358 <span class="preprocessor"></span>
<a name="l00359"></a>00359 <span class="keywordflow">return</span> res;
<a name="l00360"></a>00360 } <span class="comment">/* PolyhedronLTQ */</span>
<a name="l00361"></a>00361
<a name="l00362"></a>00362 <span class="comment">/* GaussSimplify --</span>
<a name="l00363"></a>00363 <span class="comment"> Given Mat1, a matrix of equalities, performs Gaussian elimination.</span>
<a name="l00364"></a>00364 <span class="comment"> Find a minimum basis, Returns the rank.</span>
<a name="l00365"></a>00365 <span class="comment"> Mat1 is context, Mat2 is reduced in context of Mat1</span>
<a name="l00366"></a>00366 <span class="comment">*/</span>
<a name="l00367"></a><a class="code" href="alpha_8h.html#bc2f8eafa5cd2a9e5518911f48694d50">00367</a> <span class="keywordtype">int</span> <a class="code" href="alpha_8c.html#bf25c91aa148c825f260e5c8b50fa5b1">GaussSimplify</a>(<a class="code" href="structmatrix.html">Matrix</a> *Mat1,<a class="code" href="structmatrix.html">Matrix</a> *Mat2) {
<a name="l00368"></a>00368
<a name="l00369"></a>00369 <span class="keywordtype">int</span> NbRows = Mat1-><a class="code" href="structmatrix.html#16ad614d15c6e81c0041e877b623c72d">NbRows</a>;
<a name="l00370"></a>00370 <span class="keywordtype">int</span> NbCols = Mat1-><a class="code" href="structmatrix.html#68858fd3b57684ef38bdfce13c65d182">NbColumns</a>;
<a name="l00371"></a>00371 <span class="keywordtype">int</span> *column_index;
<a name="l00372"></a>00372 <span class="keywordtype">int</span> i, j, k, <a class="code" href="polyparam_8c.html#76f11d9a0a47b94f72c2d0e77fb32240">n</a>, t, pivot, Rank;
<a name="l00373"></a>00373 Value gcd, tmp, *cp;
<a name="l00374"></a>00374
<a name="l00375"></a>00375 column_index=(<span class="keywordtype">int</span> *)malloc(NbCols * <span class="keyword">sizeof</span>(<span class="keywordtype">int</span>));
<a name="l00376"></a>00376 <span class="keywordflow">if</span> (!column_index) {
<a name="l00377"></a>00377 <a class="code" href="errormsg_8c.html#1ac398d9106568a1991a001172e0e00b">errormsg1</a>(<span class="stringliteral">"GaussSimplify"</span>, <span class="stringliteral">"outofmem"</span>, <span class="stringliteral">"out of memory space\n"</span>);
<a name="l00378"></a>00378 <a class="code" href="errormsg_8c.html#a587bf4e8c763fa5e95542df7eb070b7">Pol_status</a> = 1;
<a name="l00379"></a>00379 <span class="keywordflow">return</span> 0;
<a name="l00380"></a>00380 }
<a name="l00381"></a>00381
<a name="l00382"></a>00382 <span class="comment">/* Initialize all the 'Value' variables */</span>
<a name="l00383"></a>00383 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(gcd); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(tmp);
<a name="l00384"></a>00384
<a name="l00385"></a>00385 Rank=0;
<a name="l00386"></a>00386 <span class="keywordflow">for</span> (j=0; j<NbCols; j++) { <span class="comment">/* for each column starting at */</span>
<a name="l00387"></a>00387 <span class="keywordflow">for</span> (i=Rank; i<NbRows; i++) <span class="comment">/* diagonal, look down to find */</span>
<a name="l00388"></a>00388 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#47d32925340d2dc99ef2d4215080a60d">value_notzero_p</a>(Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i][j])) <span class="comment">/* the first non-zero entry */</span>
<a name="l00389"></a>00389 <span class="keywordflow">break</span>;
<a name="l00390"></a>00390 <span class="keywordflow">if</span> (i!=NbRows) { <span class="comment">/* was one found ? */</span>
<a name="l00391"></a>00391 <span class="keywordflow">if</span> (i!=Rank) <span class="comment">/* was it found below the diagonal?*/</span>
<a name="l00392"></a>00392 <a class="code" href="vector_8c.html#88f4d69371f16d6d28e8a87c496b0730">Vector_Exchange</a>(Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Rank],Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],NbCols);
<a name="l00393"></a>00393
<a name="l00394"></a>00394 <span class="comment">/* Normalize the pivot row */</span>
<a name="l00395"></a>00395 <a class="code" href="vector_8c.html#801486b45cadc8b9d22d99a79af91fa1">Vector_Gcd</a>(Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Rank],NbCols,&gcd);
<a name="l00396"></a>00396
<a name="l00397"></a>00397 <span class="comment">/* If (gcd >= 2) */</span>
<a name="l00398"></a>00398 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(tmp,2);
<a name="l00399"></a>00399 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#eee525b427bd616a24b5b9feb9a9f02e">value_ge</a>(gcd,tmp)) {
<a name="l00400"></a>00400 cp = Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Rank];
<a name="l00401"></a>00401 <span class="keywordflow">for</span> (k=0; k<NbCols; k++,cp++)
<a name="l00402"></a>00402 <a class="code" href="source_2arith_2arithmetique_8h.html#c633e4bcd9ed1d8412c67cda6d8ddb11">value_division</a>(*cp,*cp,gcd);
<a name="l00403"></a>00403 }
<a name="l00404"></a>00404 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#f031da52c8a160fd1ee6e4c793f8c78a">value_neg_p</a>(Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Rank][j])) {
<a name="l00405"></a>00405 cp = Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Rank];
<a name="l00406"></a>00406 <span class="keywordflow">for</span> (k=0; k<NbCols; k++,cp++)
<a name="l00407"></a>00407 <a class="code" href="source_2arith_2arithmetique_8h.html#0daf9a8ecdc14e7a274832b0d2add830">value_oppose</a>(*cp,*cp);
<a name="l00408"></a>00408 }
<a name="l00409"></a>00409 <span class="comment">/* End of normalize */</span>
<a name="l00410"></a>00410 pivot=i;
<a name="l00411"></a>00411 <span class="keywordflow">for</span> (i=0;i<NbRows;i++) <span class="comment">/* Zero out the rest of the column */</span>
<a name="l00412"></a>00412 <span class="keywordflow">if</span> (i!=Rank) {
<a name="l00413"></a>00413 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#47d32925340d2dc99ef2d4215080a60d">value_notzero_p</a>(Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i][j])) {
<a name="l00414"></a>00414 Value a, a1, a2, a1abs, a2abs;
<a name="l00415"></a>00415 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a1); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a2);
<a name="l00416"></a>00416 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a1abs); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a2abs);
<a name="l00417"></a>00417 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(a1,Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i][j]);
<a name="l00418"></a>00418 <a class="code" href="source_2arith_2arithmetique_8h.html#6eee595d027bdf64a6c529b5f138720f">value_absolute</a>(a1abs,a1);
<a name="l00419"></a>00419 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(a2,Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Rank][j]);
<a name="l00420"></a>00420 <a class="code" href="source_2arith_2arithmetique_8h.html#6eee595d027bdf64a6c529b5f138720f">value_absolute</a>(a2abs,a2);
<a name="l00421"></a>00421 <a class="code" href="source_2arith_2arithmetique_8h.html#50a34e78517e58cd1c2494a99be1e62b">value_gcd</a>(a, a1abs, a2abs);
<a name="l00422"></a>00422 <a class="code" href="source_2arith_2arithmetique_8h.html#4b01d50506e0c7691b69d57a87280f5d">value_divexact</a>(a1, a1, a);
<a name="l00423"></a>00423 <a class="code" href="source_2arith_2arithmetique_8h.html#4b01d50506e0c7691b69d57a87280f5d">value_divexact</a>(a2, a2, a);
<a name="l00424"></a>00424 <a class="code" href="source_2arith_2arithmetique_8h.html#0daf9a8ecdc14e7a274832b0d2add830">value_oppose</a>(a1,a1);
<a name="l00425"></a>00425 <a class="code" href="vector_8c.html#4700bd41bb9179b23deb7ab9d80a463d">Vector_Combine</a>(Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Rank],Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],a2,
<a name="l00426"></a>00426 a1,NbCols);
<a name="l00427"></a>00427 <a class="code" href="vector_8c.html#9e3e9df020a7c7b0ad95baffed7ed5c9">Vector_Normalize</a>(Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],NbCols);
<a name="l00428"></a>00428 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a1); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a2);
<a name="l00429"></a>00429 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a1abs); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a2abs);
<a name="l00430"></a>00430 }
<a name="l00431"></a>00431 }
<a name="l00432"></a>00432 column_index[Rank]=j;
<a name="l00433"></a>00433 Rank++;
<a name="l00434"></a>00434 }
<a name="l00435"></a>00435 } <span class="comment">/* end of Gauss elimination */</span>
<a name="l00436"></a>00436
<a name="l00437"></a>00437
<a name="l00438"></a>00438 <span class="keywordflow">if</span> (Mat2) { <span class="comment">/* Mat2 is a transformation matrix (i,j->f(i,j))....</span>
<a name="l00439"></a>00439 <span class="comment"> can't scale it because can't scale both sides of -> */</span>
<a name="l00440"></a>00440 <span class="comment">/* normalizes an affine transformation */</span>
<a name="l00441"></a>00441 <span class="comment">/* priority of forms */</span>
<a name="l00442"></a>00442 <span class="comment">/* 1. i' -> i (identity) */</span>
<a name="l00443"></a>00443 <span class="comment">/* 2. i' -> i + constant (uniform) */</span>
<a name="l00444"></a>00444 <span class="comment">/* 3. i' -> constant (broadcast) */</span>
<a name="l00445"></a>00445 <span class="comment">/* 4. i' -> j (permutation) */</span>
<a name="l00446"></a>00446 <span class="comment">/* 5. i' -> j + constant ( ) */</span>
<a name="l00447"></a>00447 <span class="comment">/* 6. i' -> i + j + constant (non-uniform) */</span>
<a name="l00448"></a>00448 <span class="keywordflow">for</span> (k=0; k<Rank; k++) {
<a name="l00449"></a>00449 j = column_index[k];
<a name="l00450"></a>00450 <span class="keywordflow">for</span> (i=0; i<(Mat2-><a class="code" href="structmatrix.html#16ad614d15c6e81c0041e877b623c72d">NbRows</a>-1);i++) { <span class="comment">/* all but the last row 0...0 1 */</span>
<a name="l00451"></a>00451 <span class="keywordflow">if</span> ((i!=j) && <a class="code" href="source_2arith_2arithmetique_8h.html#47d32925340d2dc99ef2d4215080a60d">value_notzero_p</a>(Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i][j])) {
<a name="l00452"></a>00452
<a name="l00453"></a>00453 <span class="comment">/* Remove dependency of i' on j */</span>
<a name="l00454"></a>00454 Value a, a1, a1abs, a2, a2abs;
<a name="l00455"></a>00455 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a1); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a2);
<a name="l00456"></a>00456 <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a1abs); <a class="code" href="source_2arith_2arithmetique_8h.html#f71a2ca0294a19cff0cdcbdcc052ee27">value_init</a>(a2abs);
<a name="l00457"></a>00457 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(a1,Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i][j]);
<a name="l00458"></a>00458 <a class="code" href="source_2arith_2arithmetique_8h.html#6eee595d027bdf64a6c529b5f138720f">value_absolute</a>(a1abs,a1);
<a name="l00459"></a>00459 <a class="code" href="source_2arith_2arithmetique_8h.html#864613888dc46f15679aa4f63e468f89">value_assign</a>(a2,Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[k][j]);
<a name="l00460"></a>00460 <a class="code" href="source_2arith_2arithmetique_8h.html#6eee595d027bdf64a6c529b5f138720f">value_absolute</a>(a2abs,a2);
<a name="l00461"></a>00461 <a class="code" href="source_2arith_2arithmetique_8h.html#50a34e78517e58cd1c2494a99be1e62b">value_gcd</a>(a, a1abs, a2abs);
<a name="l00462"></a>00462 <a class="code" href="source_2arith_2arithmetique_8h.html#4b01d50506e0c7691b69d57a87280f5d">value_divexact</a>(a1, a1, a);
<a name="l00463"></a>00463 <a class="code" href="source_2arith_2arithmetique_8h.html#4b01d50506e0c7691b69d57a87280f5d">value_divexact</a>(a2, a2, a);
<a name="l00464"></a>00464 <a class="code" href="source_2arith_2arithmetique_8h.html#0daf9a8ecdc14e7a274832b0d2add830">value_oppose</a>(a1,a1);
<a name="l00465"></a>00465 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#61db0ba9ba925615ec70eb92e40b9ef4">value_one_p</a>(a2)) {
<a name="l00466"></a>00466 <a class="code" href="vector_8c.html#4700bd41bb9179b23deb7ab9d80a463d">Vector_Combine</a>(Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[k],Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],a2,
<a name="l00467"></a>00467 a1,NbCols);
<a name="l00468"></a>00468
<a name="l00469"></a>00469 <span class="comment">/* Vector_Normalize(Mat2->p[i],NbCols); -- can't do T */</span>
<a name="l00470"></a>00470 } <span class="comment">/* otherwise, can't do it without mult lhs prod (2i,3j->...) */</span>
<a name="l00471"></a>00471 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a1); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a2);
<a name="l00472"></a>00472 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a1abs); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(a2abs);
<a name="l00473"></a>00473
<a name="l00474"></a>00474 }
<a name="l00475"></a>00475 <span class="keywordflow">else</span> <span class="keywordflow">if</span> ((i==j) && <a class="code" href="source_2arith_2arithmetique_8h.html#827532f2140ae2aa96e46baebae09723">value_zero_p</a>(Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i][j])) {
<a name="l00476"></a>00476
<a name="l00477"></a>00477 <span class="comment">/* 'i' does not depend on j */</span>
<a name="l00478"></a>00478 <span class="keywordflow">for</span> (n=j+1; n < (NbCols-1); n++) {
<a name="l00479"></a>00479 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#47d32925340d2dc99ef2d4215080a60d">value_notzero_p</a>(Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i][n])) { <span class="comment">/* i' depends on some n */</span>
<a name="l00480"></a>00480 <a class="code" href="source_2arith_2arithmetique_8h.html#8cc56567a4a29271559ac0fd5f6c5bfa">value_set_si</a>(tmp,1);
<a name="l00481"></a>00481 <a class="code" href="vector_8c.html#4700bd41bb9179b23deb7ab9d80a463d">Vector_Combine</a>(Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],Mat1-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[k],Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[i],tmp,
<a name="l00482"></a>00482 tmp,NbCols);
<a name="l00483"></a>00483 <span class="keywordflow">break</span>;
<a name="l00484"></a>00484 } <span class="comment">/* if 'i' depends on just a constant, then leave it alone.*/</span>
<a name="l00485"></a>00485 }
<a name="l00486"></a>00486 }
<a name="l00487"></a>00487 }
<a name="l00488"></a>00488 }
<a name="l00489"></a>00489
<a name="l00490"></a>00490 <span class="comment">/* Check last row of transformation Mat2 */</span>
<a name="l00491"></a>00491 <span class="keywordflow">for</span> (j=0; j<(NbCols-1); j++)
<a name="l00492"></a>00492 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#47d32925340d2dc99ef2d4215080a60d">value_notzero_p</a>(Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Mat2-><a class="code" href="structmatrix.html#16ad614d15c6e81c0041e877b623c72d">NbRows</a>-1][j])) {
<a name="l00493"></a>00493 <a class="code" href="errormsg_8c.html#1ac398d9106568a1991a001172e0e00b">errormsg1</a>(<span class="stringliteral">"GaussSimplify"</span>, <span class="stringliteral">"corrtrans"</span>, <span class="stringliteral">"Corrupted transformation\n"</span>);
<a name="l00494"></a>00494 <span class="keywordflow">break</span>;
<a name="l00495"></a>00495 }
<a name="l00496"></a>00496
<a name="l00497"></a>00497 <span class="keywordflow">if</span> (<a class="code" href="source_2arith_2arithmetique_8h.html#5f258141f9e07b66da01ff7f4fe3d56d">value_notone_p</a>(Mat2-><a class="code" href="structmatrix.html#2c6d840d8d911ae95c2ae4fc96f4b5ba">p</a>[Mat2-><a class="code" href="structmatrix.html#16ad614d15c6e81c0041e877b623c72d">NbRows</a>-1][NbCols-1])) {
<a name="l00498"></a>00498 <a class="code" href="errormsg_8c.html#1ac398d9106568a1991a001172e0e00b">errormsg1</a>(<span class="stringliteral">"GaussSimplify"</span>, <span class="stringliteral">"corrtrans"</span>, <span class="stringliteral">"Corrupted transformation\n"</span>);
<a name="l00499"></a>00499 }
<a name="l00500"></a>00500 }
<a name="l00501"></a>00501 <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(gcd); <a class="code" href="source_2arith_2arithmetique_8h.html#b9b282921e85a0527d462d331533d619">value_clear</a>(tmp);
<a name="l00502"></a>00502 free(column_index);
<a name="l00503"></a>00503 <span class="keywordflow">return</span> Rank;
<a name="l00504"></a>00504 } <span class="comment">/* GaussSimplify */</span>
<a name="l00505"></a>00505
<a name="l00506"></a>00506 <span class="comment">/* </span>
<a name="l00507"></a>00507 <span class="comment"> * Topologically sort 'n' polyhdera and return 0 on failure, otherwise return </span>
<a name="l00508"></a>00508 <span class="comment"> * 1 on success. Here 'L' is a an array of pointers to polyhedra, 'n' is the </span>
<a name="l00509"></a>00509 <span class="comment"> * number of polyhedra, 'index' is the level to consider for partial ordering</span>
<a name="l00510"></a>00510 <span class="comment"> * 'pdim' is the parameter space dimension, 'time' is an array of 'n' integers</span>
<a name="l00511"></a>00511 <span class="comment"> * to store logical time values, 'pvect', if not NULL, is an array of 'n' </span>
<a name="l00512"></a>00512 <span class="comment"> * integers that contains a permutation specification after call and MAXRAYS is</span>
<a name="l00513"></a>00513 <span class="comment"> * the workspace size for polyhedra operations. </span>
<a name="l00514"></a>00514 <span class="comment"> */</span>
<a name="l00515"></a><a class="code" href="alpha_8h.html#682408c0e864140d11fb7b841e5d6654">00515</a> <span class="keywordtype">int</span> <a class="code" href="alpha_8c.html#682408c0e864140d11fb7b841e5d6654">PolyhedronTSort</a> (<a class="code" href="structpolyhedron.html">Polyhedron</a> **L,<span class="keywordtype">unsigned</span> <span class="keywordtype">int</span> <a class="code" href="polyparam_8c.html#76f11d9a0a47b94f72c2d0e77fb32240">n</a>,<span class="keywordtype">unsigned</span> <span class="keywordtype">int</span> index,<span class="keywordtype">unsigned</span> <span class="keywordtype">int</span> pdim,<span class="keywordtype">int</span> *time,<span class="keywordtype">int</span> *pvect,<span class="keywordtype">unsigned</span> <span class="keywordtype">int</span> <a class="code" href="verif__ehrhart_8c.html#89fd83aa168651629c012d8655635588">MAXRAYS</a>) {
<a name="l00516"></a>00516
<a name="l00517"></a>00517 <span class="keywordtype">unsigned</span> <span class="keywordtype">int</span> <span class="keyword">const</span> nbcells = ((n*(n-1))>>1); <span class="comment">/* Number of memory cells </span>
<a name="l00518"></a>00518 <span class="comment"> to allocate, see below */</span>
<a name="l00519"></a>00519 <span class="keywordtype">int</span> *dag; <span class="comment">/* The upper triangular matrix */</span>
<a name="l00520"></a>00520 <span class="keywordtype">int</span> **<a class="code" href="vector_8c.html#a45b2e3dcf291527c5aedc420819adfc">p</a>; <span class="comment">/* Array of matrix row addresses */</span>
<a name="l00521"></a>00521 <span class="keywordtype">unsigned</span> <span class="keywordtype">int</span> i, j, k;
<a name="l00522"></a>00522 <span class="keywordtype">unsigned</span> <span class="keywordtype">int</span> t, nb, isroot;
<a name="l00523"></a>00523
<a name="l00524"></a>00524 <span class="keywordflow">if</span> (n<2) <span class="keywordflow">return</span> 0;
<a name="l00525"></a>00525
<a name="l00526"></a>00526 <span class="comment">/* we need an upper triangular matrix (example with n=4):</span>
<a name="l00527"></a>00527 <span class="comment"></span>
<a name="l00528"></a>00528 <span class="comment"> . o o o</span>
<a name="l00529"></a>00529 <span class="comment"> . . o o . are unuseful cells, o are useful cells</span>
<a name="l00530"></a>00530 <span class="comment"> . . . o</span>
<a name="l00531"></a>00531 <span class="comment"> . . . .</span>
<a name="l00532"></a>00532 <span class="comment"> </span>
<a name="l00533"></a>00533 <span class="comment"> so we need to allocate (n)(n-1)/2 cells</span>
<a name="l00534"></a>00534 <span class="comment"> - dag will point to this memory.</span>
<a name="l00535"></a>00535 <span class="comment"> - p[i] will point to row i of the matrix</span>
<a name="l00536"></a>00536 <span class="comment"> p[0] = dag - 1 (such that p[0][1] == dag[0])</span>
<a name="l00537"></a>00537 <span class="comment"> p[1] = dag - 1 + (n-1)</span>
<a name="l00538"></a>00538 <span class="comment"> p[2] = dag - 1 + (n-1) + (n-2)</span>
<a name="l00539"></a>00539 <span class="comment"> ...</span>
<a name="l00540"></a>00540 <span class="comment"> p[i] = p[i-1] + (n-1-i)</span>
<a name="l00541"></a>00541 <span class="comment"> */</span>
<a name="l00542"></a>00542
<a name="l00543"></a>00543 <span class="comment">/* malloc n(n-1)/2 for dag and n-1 for p (node n does not have any row) */</span>
<a name="l00544"></a>00544 dag = (<span class="keywordtype">int</span> *) malloc(nbcells*<span class="keyword">sizeof</span>(<span class="keywordtype">int</span>));
<a name="l00545"></a>00545 <span class="keywordflow">if</span> (!dag) <span class="keywordflow">return</span> 0;
<a name="l00546"></a>00546 p = (<span class="keywordtype">int</span> **) malloc ((n-1) * <span class="keyword">sizeof</span>(<span class="keywordtype">int</span> *));
<a name="l00547"></a>00547 <span class="keywordflow">if</span> (!p) {
<a name="l00548"></a>00548 free(dag); <span class="keywordflow">return</span> 0;
<a name="l00549"></a>00549 }
<a name="l00550"></a>00550
<a name="l00551"></a>00551 <span class="comment">/* Initialize 'p' and 'dag' */</span>
<a name="l00552"></a>00552 p[0] = dag-1;
<a name="l00553"></a>00553 <span class="keywordflow">for</span> (i=1; i<n-1; i++)
<a name="l00554"></a>00554 p[i] = p[i-1] + (n-1)-i;
<a name="l00555"></a>00555 <span class="keywordflow">for</span> (i=0; i<nbcells; i++)
<a name="l00556"></a>00556 dag[i] = -2; <span class="comment">/* -2 means 'not computed yet' */</span>
<a name="l00557"></a>00557 <span class="keywordflow">for</span> (i=0; i<n; i++) time[i] = -1;
<a name="l00558"></a>00558
<a name="l00559"></a>00559 <span class="comment">/* Compute the dag using transitivity to reduce the number of */</span>
<a name="l00560"></a>00560 <span class="comment">/* PolyhedronLTQ calls. */</span>
<a name="l00561"></a>00561 <span class="keywordflow">for</span> (i=0; i<n-1; i++) {
<a name="l00562"></a>00562 <a class="code" href="polyhedron_8h.html#729f721c450a034c1219447745c6b123">POL_ENSURE_FACETS</a>(L[i]);
<a name="l00563"></a>00563 <a class="code" href="polyhedron_8h.html#317e5eb888f5c14e1229742f31169378">POL_ENSURE_VERTICES</a>(L[i]);
<a name="l00564"></a>00564 <span class="keywordflow">for</span> (j=i+1; j<n; j++) {
<a name="l00565"></a>00565 <span class="keywordflow">if</span> (p[i][j] == -2) <span class="comment">/* not computed yes */</span>
<a name="l00566"></a>00566 p[i][j] = <a class="code" href="alpha_8c.html#0f28b3a3172da3ef981f4a57c7ead628">PolyhedronLTQ</a>(L[i], L[j], index, pdim, MAXRAYS);
<a name="l00567"></a>00567 <span class="keywordflow">if</span> (p[i][j] != 0) {
<a name="l00568"></a>00568
<a name="l00569"></a>00569 <span class="comment">/* if p[i][j] is 1 or -1, look for -p[i][j] on the same row:</span>
<a name="l00570"></a>00570 <span class="comment"> p[i][j] == -p[i][k] ==> p[j][k] = p[i][k] (transitivity)</span>
<a name="l00571"></a>00571 <span class="comment"> note: p[r][c] == -p[c][r], use this to avoid reading or writing</span>
<a name="l00572"></a>00572 <span class="comment"> under the matrix diagonal</span>
<a name="l00573"></a>00573 <span class="comment"> */</span>
<a name="l00574"></a>00574
<a name="l00575"></a>00575 <span class="comment">/* first, k<i so look for p[i][j] == p[k][i] (i.e. -p[i][k]) */</span>
<a name="l00576"></a>00576 <span class="keywordflow">for</span> (k=0; k<i; k++)
<a name="l00577"></a>00577 <span class="keywordflow">if</span> (p[k][i] == p[i][j]) p[k][j] = p[k][i];
<a name="l00578"></a>00578
<a name="l00579"></a>00579 <span class="comment">/* then, i<k<j so look for p[i][j] == -p[i][k] */</span>
<a name="l00580"></a>00580 <span class="keywordflow">for</span> (k=i+1; k<j; k++)
<a name="l00581"></a>00581 <span class="keywordflow">if</span> (p[i][k] == -p[i][j]) p[k][j] = -p[i][k];
<a name="l00582"></a>00582
<a name="l00583"></a>00583 <span class="comment">/* last, k>j same search but */</span>
<a name="l00584"></a>00584 <span class="keywordflow">for</span> (k=j+1; k<n; k++)
<a name="l00585"></a>00585 <span class="keywordflow">if</span> (p[i][k] == -p[i][j]) p[j][k] = p[i][k];
<a name="l00586"></a>00586 }
<a name="l00587"></a>00587 }
<a name="l00588"></a>00588 }
<a name="l00589"></a>00589
<a name="l00590"></a>00590 <span class="comment">/* walk thru the dag to compute the partial order (and optionally</span>
<a name="l00591"></a>00591 <span class="comment"> the permutation)</span>
<a name="l00592"></a>00592 <span class="comment"> Note: this is not the fastest way to do it but it takes</span>
<a name="l00593"></a>00593 <span class="comment"> negligible time compared to a single call of PolyhedronLTQ !</span>
<a name="l00594"></a>00594 <span class="comment"> Each macro-step (while loop) gives the same logical time to all</span>
<a name="l00595"></a>00595 <span class="comment"> found roots and optionally add these nodes in the permutation vector</span>
<a name="l00596"></a>00596 <span class="comment"> */</span>
<a name="l00597"></a>00597
<a name="l00598"></a>00598 t = 0; <span class="comment">/* current logical time, assigned to current roots and</span>
<a name="l00599"></a>00599 <span class="comment"> increased by 1 at the end of each step */</span>
<a name="l00600"></a>00600 nb = 0; <span class="comment">/* number of processed nodes (have been given a time) */</span>
<a name="l00601"></a>00601 <span class="keywordflow">while</span> (nb<n) {
<a name="l00602"></a>00602 <span class="keywordflow">for</span> (i=0; i<n; i++) { <span class="comment">/* search for roots */</span>
<a name="l00603"></a>00603 <span class="comment">/* for any node, if it's not already been given a logical time</span>
<a name="l00604"></a>00604 <span class="comment"> then walk thru the node row; if we find a 1 at some column j,</span>
<a name="l00605"></a>00605 <span class="comment"> it means that node j preceeds the current node, so it is not</span>
<a name="l00606"></a>00606 <span class="comment"> a root */</span>
<a name="l00607"></a>00607 <span class="keywordflow">if</span> (time[i]<0) {
<a name="l00608"></a>00608 isroot = 1; <span class="comment">/* assume that it is until it is definitely not */</span>
<a name="l00609"></a>00609 <span class="comment">/* first search on a column */</span>
<a name="l00610"></a>00610 <span class="keywordflow">for</span> (j=0; j<i; j++) {
<a name="l00611"></a>00611 <span class="keywordflow">if</span> (p[j][i]==-1) { <span class="comment">/* found a node lower than it */</span>
<a name="l00612"></a>00612 isroot = 0; <span class="keywordflow">break</span>;
<a name="l00613"></a>00613 }
<a name="l00614"></a>00614 }
<a name="l00615"></a>00615 <span class="keywordflow">if</span> <span class="comment">/*still*/</span> (isroot)
<a name="l00616"></a>00616 <span class="keywordflow">for</span> (j=i+1; j<n; j++) {
<a name="l00617"></a>00617 <span class="keywordflow">if</span> (p[i][j]==1) { <span class="comment">/* greater than this node */</span>
<a name="l00618"></a>00618 isroot = 0; <span class="keywordflow">break</span>;
<a name="l00619"></a>00619 }
<a name="l00620"></a>00620 }
<a name="l00621"></a>00621 <span class="keywordflow">if</span> (isroot) { <span class="comment">/* then it definitely is */</span>
<a name="l00622"></a>00622 time[i] = t; <span class="comment">/* this node gets current logical time */</span>
<a name="l00623"></a>00623 <span class="keywordflow">if</span> (pvect)
<a name="l00624"></a>00624 pvect[nb] = i+1; <span class="comment">/* nodes will be numbered from 1 to n */</span>
<a name="l00625"></a>00625 nb++; <span class="comment">/* one more node processed */</span>
<a name="l00626"></a>00626 }
<a name="l00627"></a>00627 }
<a name="l00628"></a>00628 }
<a name="l00629"></a>00629 <span class="comment">/* now make roots become neutral, i.e. non comparable with all other nodes */</span>
<a name="l00630"></a>00630 <span class="keywordflow">for</span> (i=0; i<n; i++) {
<a name="l00631"></a>00631 <span class="keywordflow">if</span> (time[i]==t) {
<a name="l00632"></a>00632 <span class="keywordflow">for</span> (j=0; j<i; j++)
<a name="l00633"></a>00633 p[j][i] = 0;
<a name="l00634"></a>00634 <span class="keywordflow">for</span> (j=i+1; j<n; j++)
<a name="l00635"></a>00635 p[i][j] = 0;
<a name="l00636"></a>00636 }
<a name="l00637"></a>00637 }
<a name="l00638"></a>00638 t++; <span class="comment">/* ready for next set of root nodes */</span>
<a name="l00639"></a>00639 }
<a name="l00640"></a>00640
<a name="l00641"></a>00641 free (p); <span class="comment">/* let's be clean, collect the garbage */</span>
<a name="l00642"></a>00642 free (dag);
<a name="l00643"></a>00643 <span class="keywordflow">return</span> 1;
<a name="l00644"></a>00644 } <span class="comment">/* PolyhedronTSort */</span>
<a name="l00645"></a>00645
<a name="l00646"></a>00646
<a name="l00647"></a>00647
<a name="l00648"></a>00648
</pre></div></div>
<hr size="1"><address style="text-align: right;"><small>Generated on Tue Sep 15 18:33:58 2009 for polylib by
<a href="http://www.doxygen.org/index.html">
<img src="doxygen.png" alt="doxygen" align="middle" border="0"></a> 1.5.6 </small></address>
</body>
</html>
|
