add isl_union_set_compute_schedule
[platform/upstream/isl.git] / isl_schedule.c
1 /*
2  * Copyright 2011      INRIA Saclay
3  *
4  * Use of this software is governed by the GNU LGPLv2.1 license
5  *
6  * Written by Sven Verdoolaege, INRIA Saclay - Ile-de-France,
7  * Parc Club Orsay Universite, ZAC des vignes, 4 rue Jacques Monod,
8  * 91893 Orsay, France
9  */
10
11 #include <isl_ctx_private.h>
12 #include <isl_map_private.h>
13 #include <isl_dim_private.h>
14 #include <isl/hash.h>
15 #include <isl/constraint.h>
16 #include <isl/schedule.h>
17 #include <isl_mat_private.h>
18 #include <isl/set.h>
19 #include <isl/seq.h>
20 #include <isl_tab.h>
21 #include <isl_dim_map.h>
22 #include <isl_hmap_map_basic_set.h>
23 #include <isl_qsort.h>
24
25 /*
26  * The scheduling algorithm implemented in this file was inspired by
27  * Bondhugula et al., "Automatic Transformations for Communication-Minimized
28  * Parallelization and Locality Optimization in the Polyhedral Model".
29  */
30
31
32 /* The schedule for an individual domain, plus information about the bands.
33  * In particular, we keep track of the number of bands and for each
34  * band, the starting position of the next band.  The first band starts at
35  * position 0.
36  */
37 struct isl_schedule_node {
38         isl_map *sched;
39         int      n_band;
40         int     *band_end;
41 };
42
43 /* Information about the computed schedule.
44  * n is the number of nodes/domains/statements.
45  * n_band is the maximal number of bands.
46  * n_total_row is the number of coordinates of the schedule.
47  * dim contains a description of the parameters.
48  */
49 struct isl_schedule {
50         int n;
51         int n_band;
52         int n_total_row;
53         isl_dim *dim;
54
55         struct isl_schedule_node node[1];
56 };
57
58 /* Internal information about a node that is used during the construction
59  * of a schedule.
60  * dim represents the space in which the domain lives
61  * sched is a matrix representation of the schedule being constructed
62  *      for this node
63  * sched_map is an isl_map representation of the same (partial) schedule
64  *      sched_map may be NULL
65  * rank is the number of linearly independent rows in the linear part
66  *      of sched
67  * the columns of cmap represent a change of basis for the schedule
68  *      coefficients; the first rank columns span the linear part of
69  *      the schedule rows
70  * start is the first variable in the LP problem in the sequences that
71  *      represents the schedule coefficients of this node
72  * nvar is the dimension of the domain
73  * nparam is the number of parameters or 0 if we are not constructing
74  *      a parametric schedule
75  *
76  * scc is the index of SCC (or WCC) this node belongs to
77  *
78  * band contains the band index for each of the rows of the schedule
79  *
80  * index, min_index and on_stack are used during the SCC detection
81  * index represents the order in which nodes are visited.
82  * min_index is the index of the root of a (sub)component.
83  * on_stack indicates whether the node is currently on the stack.
84  */
85 struct isl_sched_node {
86         isl_dim *dim;
87         isl_mat *sched;
88         isl_map *sched_map;
89         int      rank;
90         isl_mat *cmap;
91         int      start;
92         int      nvar;
93         int      nparam;
94
95         int      scc;
96
97         int     *band;
98
99         /* scc detection */
100         int      index;
101         int      min_index;
102         int      on_stack;
103 };
104
105 static int node_has_dim(const void *entry, const void *val)
106 {
107         struct isl_sched_node *node = (struct isl_sched_node *)entry;
108         isl_dim *dim = (isl_dim *)val;
109
110         return isl_dim_equal(node->dim, dim);
111 }
112
113 /* An edge in the dependence graph.  An edge may be used to
114  * ensure validity of the generated schedule, to minimize the dependence
115  * distance or both
116  *
117  * map is the dependence relation
118  * src is the source node
119  * dst is the sink node
120  * validity is set if the edge is used to ensure correctness
121  * proximity is set if the edge is used to minimize dependence distances
122  *
123  * For validity edges, start and end mark the sequence of inequality
124  * constraints in the LP problem that encode the validity constraint
125  * corresponding to this edge.
126  */
127 struct isl_sched_edge {
128         isl_map *map;
129
130         struct isl_sched_node *src;
131         struct isl_sched_node *dst;
132
133         int validity;
134         int proximity;
135
136         int start;
137         int end;
138 };
139
140 /* Internal information about the dependence graph used during
141  * the construction of the schedule.
142  *
143  * intra_hmap is a cache, mapping dependence relations to their dual,
144  *      for dependences from a node to itself
145  * inter_hmap is a cache, mapping dependence relations to their dual,
146  *      for dependences between distinct nodes
147  *
148  * n is the number of nodes
149  * node is the list of nodes
150  * maxvar is the maximal number of variables over all nodes
151  * n_row is the current (maximal) number of linearly independent
152  *      rows in the node schedules
153  * n_total_row is the current number of rows in the node schedules
154  * n_band is the current number of completed bands
155  * band_start is the starting row in the node schedules of the current band
156  * root is set if this graph is the original dependence graph,
157  *      without any splitting
158  *
159  * sorted contains a list of node indices sorted according to the
160  *      SCC to which a node belongs
161  *
162  * n_edge is the number of edges
163  * edge is the list of edges
164  * edge_table contains pointers into the edge array, hashed on the source
165  *      and sink spaces; the table only contains edges that represent
166  *      validity constraints (and that may or may not also represent proximity
167  *      constraints)
168  *
169  * node_table contains pointers into the node array, hashed on the space
170  *
171  * region contains a list of variable sequences that should be non-trivial
172  *
173  * lp contains the (I)LP problem used to obtain new schedule rows
174  *
175  * src_scc and dst_scc are the source and sink SCCs of an edge with
176  *      conflicting constraints
177  *
178  * scc, sp, index and stack are used during the detection of SCCs
179  * scc is the number of the next SCC
180  * stack contains the nodes on the path from the root to the current node
181  * sp is the stack pointer
182  * index is the index of the last node visited
183  */
184 struct isl_sched_graph {
185         isl_hmap_map_basic_set *intra_hmap;
186         isl_hmap_map_basic_set *inter_hmap;
187
188         struct isl_sched_node *node;
189         int n;
190         int maxvar;
191         int n_row;
192
193         int *sorted;
194
195         int n_band;
196         int n_total_row;
197         int band_start;
198
199         int root;
200
201         struct isl_sched_edge *edge;
202         int n_edge;
203         struct isl_hash_table *edge_table;
204
205         struct isl_hash_table *node_table;
206         struct isl_region *region;
207
208         isl_basic_set *lp;
209
210         int src_scc;
211         int dst_scc;
212
213         /* scc detection */
214         int scc;
215         int sp;
216         int index;
217         int *stack;
218 };
219
220 /* Initialize node_table based on the list of nodes.
221  */
222 static int graph_init_table(isl_ctx *ctx, struct isl_sched_graph *graph)
223 {
224         int i;
225
226         graph->node_table = isl_hash_table_alloc(ctx, graph->n);
227         if (!graph->node_table)
228                 return -1;
229
230         for (i = 0; i < graph->n; ++i) {
231                 struct isl_hash_table_entry *entry;
232                 uint32_t hash;
233
234                 hash = isl_dim_get_hash(graph->node[i].dim);
235                 entry = isl_hash_table_find(ctx, graph->node_table, hash,
236                                             &node_has_dim,
237                                             graph->node[i].dim, 1);
238                 if (!entry)
239                         return -1;
240                 entry->data = &graph->node[i];
241         }
242
243         return 0;
244 }
245
246 /* Return a pointer to the node that lives within the given space,
247  * or NULL if there is no such node.
248  */
249 static struct isl_sched_node *graph_find_node(isl_ctx *ctx,
250         struct isl_sched_graph *graph, __isl_keep isl_dim *dim)
251 {
252         struct isl_hash_table_entry *entry;
253         uint32_t hash;
254
255         hash = isl_dim_get_hash(dim);
256         entry = isl_hash_table_find(ctx, graph->node_table, hash,
257                                     &node_has_dim, dim, 0);
258
259         return entry ? entry->data : NULL;
260 }
261
262 static int edge_has_src_and_dst(const void *entry, const void *val)
263 {
264         const struct isl_sched_edge *edge = entry;
265         const struct isl_sched_edge *temp = val;
266
267         return edge->src == temp->src && edge->dst == temp->dst;
268 }
269
270 /* Initialize edge_table based on the list of edges.
271  * Only edges with validity set are added to the table.
272  */
273 static int graph_init_edge_table(isl_ctx *ctx, struct isl_sched_graph *graph)
274 {
275         int i;
276
277         graph->edge_table = isl_hash_table_alloc(ctx, graph->n_edge);
278         if (!graph->edge_table)
279                 return -1;
280
281         for (i = 0; i < graph->n_edge; ++i) {
282                 struct isl_hash_table_entry *entry;
283                 uint32_t hash;
284
285                 if (!graph->edge[i].validity)
286                         continue;
287
288                 hash = isl_hash_init();
289                 hash = isl_hash_builtin(hash, graph->edge[i].src);
290                 hash = isl_hash_builtin(hash, graph->edge[i].dst);
291                 entry = isl_hash_table_find(ctx, graph->edge_table, hash,
292                                             &edge_has_src_and_dst,
293                                             &graph->edge[i], 1);
294                 if (!entry)
295                         return -1;
296                 entry->data = &graph->edge[i];
297         }
298
299         return 0;
300 }
301
302 /* Check whether the dependence graph has a (validity) edge
303  * between the given two nodes.
304  */
305 static int graph_has_edge(struct isl_sched_graph *graph,
306         struct isl_sched_node *src, struct isl_sched_node *dst)
307 {
308         isl_ctx *ctx = isl_dim_get_ctx(src->dim);
309         struct isl_hash_table_entry *entry;
310         uint32_t hash;
311         struct isl_sched_edge temp = { .src = src, .dst = dst };
312         struct isl_sched_edge *edge;
313         int empty;
314
315         hash = isl_hash_init();
316         hash = isl_hash_builtin(hash, temp.src);
317         hash = isl_hash_builtin(hash, temp.dst);
318         entry = isl_hash_table_find(ctx, graph->edge_table, hash,
319                                     &edge_has_src_and_dst, &temp, 0);
320         if (!entry)
321                 return 0;
322
323         edge = entry->data;
324         empty = isl_map_fast_is_empty(edge->map);
325         if (empty < 0)
326                 return -1;
327
328         return !empty;
329 }
330
331 static int graph_alloc(isl_ctx *ctx, struct isl_sched_graph *graph,
332         int n_node, int n_edge)
333 {
334         int i;
335
336         graph->n = n_node;
337         graph->n_edge = n_edge;
338         graph->node = isl_calloc_array(ctx, struct isl_sched_node, graph->n);
339         graph->sorted = isl_calloc_array(ctx, int, graph->n);
340         graph->region = isl_alloc_array(ctx, struct isl_region, graph->n);
341         graph->stack = isl_alloc_array(ctx, int, graph->n);
342         graph->edge = isl_calloc_array(ctx,
343                                         struct isl_sched_edge, graph->n_edge);
344
345         graph->intra_hmap = isl_hmap_map_basic_set_alloc(ctx, 2 * n_edge);
346         graph->inter_hmap = isl_hmap_map_basic_set_alloc(ctx, 2 * n_edge);
347
348         if (!graph->node || !graph->region || !graph->stack || !graph->edge ||
349             !graph->sorted)
350                 return -1;
351
352         for(i = 0; i < graph->n; ++i)
353                 graph->sorted[i] = i;
354
355         return 0;
356 }
357
358 static void graph_free(isl_ctx *ctx, struct isl_sched_graph *graph)
359 {
360         int i;
361
362         isl_hmap_map_basic_set_free(ctx, graph->intra_hmap);
363         isl_hmap_map_basic_set_free(ctx, graph->inter_hmap);
364
365         for (i = 0; i < graph->n; ++i) {
366                 isl_dim_free(graph->node[i].dim);
367                 isl_mat_free(graph->node[i].sched);
368                 isl_map_free(graph->node[i].sched_map);
369                 isl_mat_free(graph->node[i].cmap);
370                 if (graph->root)
371                         free(graph->node[i].band);
372         }
373         free(graph->node);
374         free(graph->sorted);
375         for (i = 0; i < graph->n_edge; ++i)
376                 isl_map_free(graph->edge[i].map);
377         free(graph->edge);
378         free(graph->region);
379         free(graph->stack);
380         isl_hash_table_free(ctx, graph->edge_table);
381         isl_hash_table_free(ctx, graph->node_table);
382         isl_basic_set_free(graph->lp);
383 }
384
385 /* Add a new node to the graph representing the given set.
386  */
387 static int extract_node(__isl_take isl_set *set, void *user)
388 {
389         int i;
390         int nvar, nparam;
391         isl_ctx *ctx;
392         isl_dim *dim;
393         isl_mat *sched;
394         struct isl_sched_graph *graph = user;
395         int *band;
396
397         ctx = isl_set_get_ctx(set);
398         dim = isl_set_get_dim(set);
399         isl_set_free(set);
400         nvar = isl_dim_size(dim, isl_dim_set);
401         nparam = isl_dim_size(dim, isl_dim_param);
402         if (!ctx->opt->schedule_parametric)
403                 nparam = 0;
404         sched = isl_mat_alloc(ctx, 0, 1 + nparam + nvar);
405         graph->node[graph->n].dim = dim;
406         graph->node[graph->n].nvar = nvar;
407         graph->node[graph->n].nparam = nparam;
408         graph->node[graph->n].sched = sched;
409         graph->node[graph->n].sched_map = NULL;
410         band = isl_alloc_array(ctx, int, graph->n_edge + nvar);
411         graph->node[graph->n].band = band;
412         graph->n++;
413
414         if (!sched || !band)
415                 return -1;
416
417         return 0;
418 }
419
420 /* Add a new edge to the graph based on the given map.
421  * Edges are first extracted from the validity dependences,
422  * from which the edge_table is constructed.
423  * Afterwards, the proximity dependences are added.  If a proximity
424  * dependence relation happens to be identical to one of the
425  * validity dependence relations added before, then we don't create
426  * a new edge, but instead mark the original edge as also representing
427  * a proximity dependence.
428  */
429 static int extract_edge(__isl_take isl_map *map, void *user)
430 {
431         isl_ctx *ctx = isl_map_get_ctx(map);
432         struct isl_sched_graph *graph = user;
433         struct isl_sched_node *src, *dst;
434         isl_dim *dim;
435
436         dim = isl_dim_domain(isl_map_get_dim(map));
437         src = graph_find_node(ctx, graph, dim);
438         isl_dim_free(dim);
439         dim = isl_dim_range(isl_map_get_dim(map));
440         dst = graph_find_node(ctx, graph, dim);
441         isl_dim_free(dim);
442
443         if (!src || !dst) {
444                 isl_map_free(map);
445                 return 0;
446         }
447
448         graph->edge[graph->n_edge].src = src;
449         graph->edge[graph->n_edge].dst = dst;
450         graph->edge[graph->n_edge].map = map;
451         graph->edge[graph->n_edge].validity = !graph->edge_table;
452         graph->edge[graph->n_edge].proximity = !!graph->edge_table;
453         graph->n_edge++;
454
455         if (graph->edge_table) {
456                 uint32_t hash;
457                 struct isl_hash_table_entry *entry;
458                 struct isl_sched_edge *edge;
459                 int is_equal;
460
461                 hash = isl_hash_init();
462                 hash = isl_hash_builtin(hash, src);
463                 hash = isl_hash_builtin(hash, dst);
464                 entry = isl_hash_table_find(ctx, graph->edge_table, hash,
465                                             &edge_has_src_and_dst,
466                                             &graph->edge[graph->n_edge - 1], 0);
467                 if (!entry)
468                         return 0;
469                 edge = entry->data;
470                 is_equal = isl_map_fast_is_equal(map, edge->map);
471                 if (is_equal < 0)
472                         return -1;
473                 if (!is_equal)
474                         return 0;
475
476                 graph->n_edge--;
477                 edge->proximity = 1;
478                 isl_map_free(map);
479         }
480
481         return 0;
482 }
483
484 /* Check whether there is a validity dependence from src to dst,
485  * forcing dst to follow src.
486  */
487 static int node_follows(struct isl_sched_graph *graph, 
488         struct isl_sched_node *dst, struct isl_sched_node *src)
489 {
490         return graph_has_edge(graph, src, dst);
491 }
492
493 /* Perform Tarjan's algorithm for computing the strongly connected components
494  * in the dependence graph (only validity edges).
495  * If directed is not set, we consider the graph to be undirected and
496  * we effectively compute the (weakly) connected components.
497  */
498 static int detect_sccs_tarjan(struct isl_sched_graph *g, int i, int directed)
499 {
500         int j;
501
502         g->node[i].index = g->index;
503         g->node[i].min_index = g->index;
504         g->node[i].on_stack = 1;
505         g->index++;
506         g->stack[g->sp++] = i;
507
508         for (j = g->n - 1; j >= 0; --j) {
509                 int f;
510
511                 if (j == i)
512                         continue;
513                 if (g->node[j].index >= 0 &&
514                         (!g->node[j].on_stack ||
515                          g->node[j].index > g->node[i].min_index))
516                         continue;
517                 
518                 f = node_follows(g, &g->node[i], &g->node[j]);
519                 if (f < 0)
520                         return -1;
521                 if (!f && !directed) {
522                         f = node_follows(g, &g->node[j], &g->node[i]);
523                         if (f < 0)
524                                 return -1;
525                 }
526                 if (!f)
527                         continue;
528                 if (g->node[j].index < 0) {
529                         detect_sccs_tarjan(g, j, directed);
530                         if (g->node[j].min_index < g->node[i].min_index)
531                                 g->node[i].min_index = g->node[j].min_index;
532                 } else if (g->node[j].index < g->node[i].min_index)
533                         g->node[i].min_index = g->node[j].index;
534         }
535
536         if (g->node[i].index != g->node[i].min_index)
537                 return 0;
538
539         do {
540                 j = g->stack[--g->sp];
541                 g->node[j].on_stack = 0;
542                 g->node[j].scc = g->scc;
543         } while (j != i);
544         g->scc++;
545
546         return 0;
547 }
548
549 static int detect_ccs(struct isl_sched_graph *graph, int directed)
550 {
551         int i;
552
553         graph->index = 0;
554         graph->sp = 0;
555         graph->scc = 0;
556         for (i = graph->n - 1; i >= 0; --i)
557                 graph->node[i].index = -1;
558
559         for (i = graph->n - 1; i >= 0; --i) {
560                 if (graph->node[i].index >= 0)
561                         continue;
562                 if (detect_sccs_tarjan(graph, i, directed) < 0)
563                         return -1;
564         }
565
566         return 0;
567 }
568
569 /* Apply Tarjan's algorithm to detect the strongly connected components
570  * in the dependence graph.
571  */
572 static int detect_sccs(struct isl_sched_graph *graph)
573 {
574         return detect_ccs(graph, 1);
575 }
576
577 /* Apply Tarjan's algorithm to detect the (weakly) connected components
578  * in the dependence graph.
579  */
580 static int detect_wccs(struct isl_sched_graph *graph)
581 {
582         return detect_ccs(graph, 0);
583 }
584
585 static int cmp_scc(const void *a, const void *b, void *data)
586 {
587         struct isl_sched_graph *graph = data;
588         const int *i1 = a;
589         const int *i2 = b;
590
591         return graph->node[*i1].scc - graph->node[*i2].scc;
592 }
593
594 /* Sort the elements of graph->sorted according to the corresponding SCCs.
595  */
596 static void sort_sccs(struct isl_sched_graph *graph)
597 {
598         isl_quicksort(graph->sorted, graph->n, sizeof(int), &cmp_scc, graph);
599 }
600
601 /* Given a dependence relation R from a node to itself,
602  * construct the set of coefficients of valid constraints for elements
603  * in that dependence relation.
604  * In particular, the result contains tuples of coefficients
605  * c_0, c_n, c_x such that
606  *
607  *      c_0 + c_n n + c_x y - c_x x >= 0 for each (x,y) in R
608  *
609  * or, equivalently,
610  *
611  *      c_0 + c_n n + c_x d >= 0 for each d in delta R = { y - x | (x,y) in R }
612  *
613  * We choose here to compute the dual of delta R.
614  * Alternatively, we could have computed the dual of R, resulting
615  * in a set of tuples c_0, c_n, c_x, c_y, and then
616  * plugged in (c_0, c_n, c_x, -c_x).
617  */
618 static __isl_give isl_basic_set *intra_coefficients(
619         struct isl_sched_graph *graph, __isl_take isl_map *map)
620 {
621         isl_ctx *ctx = isl_map_get_ctx(map);
622         isl_set *delta;
623         isl_basic_set *coef;
624
625         if (isl_hmap_map_basic_set_has(ctx, graph->intra_hmap, map))
626                 return isl_hmap_map_basic_set_get(ctx, graph->intra_hmap, map);
627
628         delta = isl_set_remove_divs(isl_map_deltas(isl_map_copy(map)));
629         coef = isl_set_coefficients(delta);
630         isl_hmap_map_basic_set_set(ctx, graph->intra_hmap, map,
631                                         isl_basic_set_copy(coef));
632
633         return coef;
634 }
635
636 /* Given a dependence relation R, * construct the set of coefficients
637  * of valid constraints for elements in that dependence relation.
638  * In particular, the result contains tuples of coefficients
639  * c_0, c_n, c_x, c_y such that
640  *
641  *      c_0 + c_n n + c_x x + c_y y >= 0 for each (x,y) in R
642  *
643  */
644 static __isl_give isl_basic_set *inter_coefficients(
645         struct isl_sched_graph *graph, __isl_take isl_map *map)
646 {
647         isl_ctx *ctx = isl_map_get_ctx(map);
648         isl_set *set;
649         isl_basic_set *coef;
650
651         if (isl_hmap_map_basic_set_has(ctx, graph->inter_hmap, map))
652                 return isl_hmap_map_basic_set_get(ctx, graph->inter_hmap, map);
653
654         set = isl_map_wrap(isl_map_remove_divs(isl_map_copy(map)));
655         coef = isl_set_coefficients(set);
656         isl_hmap_map_basic_set_set(ctx, graph->inter_hmap, map,
657                                         isl_basic_set_copy(coef));
658
659         return coef;
660 }
661
662 /* Add constraints to graph->lp that force validity for the given
663  * dependence from a node i to itself.
664  * That is, add constraints that enforce
665  *
666  *      (c_i_0 + c_i_n n + c_i_x y) - (c_i_0 + c_i_n n + c_i_x x)
667  *      = c_i_x (y - x) >= 0
668  *
669  * for each (x,y) in R.
670  * We obtain general constraints on coefficients (c_0, c_n, c_x)
671  * of valid constraints for (y - x) and then plug in (0, 0, c_i_x^+ - c_i_x^-),
672  * where c_i_x = c_i_x^+ - c_i_x^-, with c_i_x^+ and c_i_x^- non-negative.
673  * In graph->lp, the c_i_x^- appear before their c_i_x^+ counterpart.
674  *
675  * Actually, we do not construct constraints for the c_i_x themselves,
676  * but for the coefficients of c_i_x written as a linear combination
677  * of the columns in node->cmap.
678  */
679 static int add_intra_validity_constraints(struct isl_sched_graph *graph,
680         struct isl_sched_edge *edge)
681 {
682         unsigned total;
683         isl_map *map = isl_map_copy(edge->map);
684         isl_ctx *ctx = isl_map_get_ctx(map);
685         isl_dim *dim;
686         isl_dim_map *dim_map;
687         isl_basic_set *coef;
688         struct isl_sched_node *node = edge->src;
689
690         coef = intra_coefficients(graph, map);
691
692         dim = isl_dim_domain(isl_dim_unwrap(isl_basic_set_get_dim(coef)));
693
694         coef = isl_basic_set_transform_dims(coef, isl_dim_set,
695                     isl_dim_size(dim, isl_dim_set), isl_mat_copy(node->cmap));
696
697         total = isl_basic_set_total_dim(graph->lp);
698         dim_map = isl_dim_map_alloc(ctx, total);
699         isl_dim_map_range(dim_map, node->start + 2 * node->nparam + 1, 2,
700                           isl_dim_size(dim, isl_dim_set), 1,
701                           node->nvar, -1);
702         isl_dim_map_range(dim_map, node->start + 2 * node->nparam + 2, 2,
703                           isl_dim_size(dim, isl_dim_set), 1,
704                           node->nvar, 1);
705         graph->lp = isl_basic_set_extend_constraints(graph->lp,
706                         coef->n_eq, coef->n_ineq);
707         graph->lp = isl_basic_set_add_constraints_dim_map(graph->lp,
708                                                            coef, dim_map);
709         isl_dim_free(dim);
710
711         return 0;
712 }
713
714 /* Add constraints to graph->lp that force validity for the given
715  * dependence from node i to node j.
716  * That is, add constraints that enforce
717  *
718  *      (c_j_0 + c_j_n n + c_j_x y) - (c_i_0 + c_i_n n + c_i_x x) >= 0
719  *
720  * for each (x,y) in R.
721  * We obtain general constraints on coefficients (c_0, c_n, c_x, c_y)
722  * of valid constraints for R and then plug in
723  * (c_j_0 - c_i_0, c_j_n^+ - c_j_n^- - (c_i_n^+ - c_i_n^-),
724  *  c_j_x^+ - c_j_x^- - (c_i_x^+ - c_i_x^-)),
725  * where c_* = c_*^+ - c_*^-, with c_*^+ and c_*^- non-negative.
726  * In graph->lp, the c_*^- appear before their c_*^+ counterpart.
727  *
728  * Actually, we do not construct constraints for the c_*_x themselves,
729  * but for the coefficients of c_*_x written as a linear combination
730  * of the columns in node->cmap.
731  */
732 static int add_inter_validity_constraints(struct isl_sched_graph *graph,
733         struct isl_sched_edge *edge)
734 {
735         unsigned total;
736         isl_map *map = isl_map_copy(edge->map);
737         isl_ctx *ctx = isl_map_get_ctx(map);
738         isl_dim *dim;
739         isl_dim_map *dim_map;
740         isl_basic_set *coef;
741         struct isl_sched_node *src = edge->src;
742         struct isl_sched_node *dst = edge->dst;
743
744         coef = inter_coefficients(graph, map);
745
746         dim = isl_dim_domain(isl_dim_unwrap(isl_basic_set_get_dim(coef)));
747
748         coef = isl_basic_set_transform_dims(coef, isl_dim_set,
749                     isl_dim_size(dim, isl_dim_set), isl_mat_copy(src->cmap));
750         coef = isl_basic_set_transform_dims(coef, isl_dim_set,
751                     isl_dim_size(dim, isl_dim_set) + src->nvar,
752                     isl_mat_copy(dst->cmap));
753
754         total = isl_basic_set_total_dim(graph->lp);
755         dim_map = isl_dim_map_alloc(ctx, total);
756
757         isl_dim_map_range(dim_map, dst->start, 0, 0, 0, 1, 1);
758         isl_dim_map_range(dim_map, dst->start + 1, 2, 1, 1, dst->nparam, -1);
759         isl_dim_map_range(dim_map, dst->start + 2, 2, 1, 1, dst->nparam, 1);
760         isl_dim_map_range(dim_map, dst->start + 2 * dst->nparam + 1, 2,
761                           isl_dim_size(dim, isl_dim_set) + src->nvar, 1,
762                           dst->nvar, -1);
763         isl_dim_map_range(dim_map, dst->start + 2 * dst->nparam + 2, 2,
764                           isl_dim_size(dim, isl_dim_set) + src->nvar, 1,
765                           dst->nvar, 1);
766
767         isl_dim_map_range(dim_map, src->start, 0, 0, 0, 1, -1);
768         isl_dim_map_range(dim_map, src->start + 1, 2, 1, 1, src->nparam, 1);
769         isl_dim_map_range(dim_map, src->start + 2, 2, 1, 1, src->nparam, -1);
770         isl_dim_map_range(dim_map, src->start + 2 * src->nparam + 1, 2,
771                           isl_dim_size(dim, isl_dim_set), 1,
772                           src->nvar, 1);
773         isl_dim_map_range(dim_map, src->start + 2 * src->nparam + 2, 2,
774                           isl_dim_size(dim, isl_dim_set), 1,
775                           src->nvar, -1);
776
777         edge->start = graph->lp->n_ineq;
778         graph->lp = isl_basic_set_extend_constraints(graph->lp,
779                         coef->n_eq, coef->n_ineq);
780         graph->lp = isl_basic_set_add_constraints_dim_map(graph->lp,
781                                                            coef, dim_map);
782         isl_dim_free(dim);
783         edge->end = graph->lp->n_ineq;
784
785         return 0;
786 }
787
788 /* Add constraints to graph->lp that bound the dependence distance for the given
789  * dependence from a node i to itself.
790  * If s = 1, we add the constraint
791  *
792  *      c_i_x (y - x) <= m_0 + m_n n
793  *
794  * or
795  *
796  *      -c_i_x (y - x) + m_0 + m_n n >= 0
797  *
798  * for each (x,y) in R.
799  * If s = -1, we add the constraint
800  *
801  *      -c_i_x (y - x) <= m_0 + m_n n
802  *
803  * or
804  *
805  *      c_i_x (y - x) + m_0 + m_n n >= 0
806  *
807  * for each (x,y) in R.
808  * We obtain general constraints on coefficients (c_0, c_n, c_x)
809  * of valid constraints for (y - x) and then plug in (m_0, m_n, -s * c_i_x),
810  * with each coefficient (except m_0) represented as a pair of non-negative
811  * coefficients.
812  *
813  * Actually, we do not construct constraints for the c_i_x themselves,
814  * but for the coefficients of c_i_x written as a linear combination
815  * of the columns in node->cmap.
816  */
817 static int add_intra_proximity_constraints(struct isl_sched_graph *graph,
818         struct isl_sched_edge *edge, int s)
819 {
820         unsigned total;
821         unsigned nparam;
822         isl_map *map = isl_map_copy(edge->map);
823         isl_ctx *ctx = isl_map_get_ctx(map);
824         isl_dim *dim;
825         isl_dim_map *dim_map;
826         isl_set *delta;
827         isl_basic_set *coef;
828         struct isl_sched_node *node = edge->src;
829
830         coef = intra_coefficients(graph, map);
831
832         dim = isl_dim_domain(isl_dim_unwrap(isl_basic_set_get_dim(coef)));
833
834         coef = isl_basic_set_transform_dims(coef, isl_dim_set,
835                     isl_dim_size(dim, isl_dim_set), isl_mat_copy(node->cmap));
836
837         nparam = isl_dim_size(node->dim, isl_dim_param);
838         total = isl_basic_set_total_dim(graph->lp);
839         dim_map = isl_dim_map_alloc(ctx, total);
840         isl_dim_map_range(dim_map, 1, 0, 0, 0, 1, 1);
841         isl_dim_map_range(dim_map, 4, 2, 1, 1, nparam, -1);
842         isl_dim_map_range(dim_map, 5, 2, 1, 1, nparam, 1);
843         isl_dim_map_range(dim_map, node->start + 2 * node->nparam + 1, 2,
844                           isl_dim_size(dim, isl_dim_set), 1,
845                           node->nvar, s);
846         isl_dim_map_range(dim_map, node->start + 2 * node->nparam + 2, 2,
847                           isl_dim_size(dim, isl_dim_set), 1,
848                           node->nvar, -s);
849         graph->lp = isl_basic_set_extend_constraints(graph->lp,
850                         coef->n_eq, coef->n_ineq);
851         graph->lp = isl_basic_set_add_constraints_dim_map(graph->lp,
852                                                            coef, dim_map);
853         isl_dim_free(dim);
854
855         return 0;
856 }
857
858 /* Add constraints to graph->lp that bound the dependence distance for the given
859  * dependence from node i to node j.
860  * If s = 1, we add the constraint
861  *
862  *      (c_j_0 + c_j_n n + c_j_x y) - (c_i_0 + c_i_n n + c_i_x x)
863  *              <= m_0 + m_n n
864  *
865  * or
866  *
867  *      -(c_j_0 + c_j_n n + c_j_x y) + (c_i_0 + c_i_n n + c_i_x x) +
868  *              m_0 + m_n n >= 0
869  *
870  * for each (x,y) in R.
871  * If s = -1, we add the constraint
872  *
873  *      -((c_j_0 + c_j_n n + c_j_x y) - (c_i_0 + c_i_n n + c_i_x x))
874  *              <= m_0 + m_n n
875  *
876  * or
877  *
878  *      (c_j_0 + c_j_n n + c_j_x y) - (c_i_0 + c_i_n n + c_i_x x) +
879  *              m_0 + m_n n >= 0
880  *
881  * for each (x,y) in R.
882  * We obtain general constraints on coefficients (c_0, c_n, c_x, c_y)
883  * of valid constraints for R and then plug in
884  * (m_0 - s*c_j_0 + s*c_i_0, m_n - s*c_j_n + s*c_i_n,
885  *  -s*c_j_x+s*c_i_x)
886  * with each coefficient (except m_0, c_j_0 and c_i_0)
887  * represented as a pair of non-negative coefficients.
888  *
889  * Actually, we do not construct constraints for the c_*_x themselves,
890  * but for the coefficients of c_*_x written as a linear combination
891  * of the columns in node->cmap.
892  */
893 static int add_inter_proximity_constraints(struct isl_sched_graph *graph,
894         struct isl_sched_edge *edge, int s)
895 {
896         unsigned total;
897         unsigned nparam;
898         isl_map *map = isl_map_copy(edge->map);
899         isl_ctx *ctx = isl_map_get_ctx(map);
900         isl_dim *dim;
901         isl_dim_map *dim_map;
902         isl_basic_set *coef;
903         struct isl_sched_node *src = edge->src;
904         struct isl_sched_node *dst = edge->dst;
905
906         coef = inter_coefficients(graph, map);
907
908         dim = isl_dim_domain(isl_dim_unwrap(isl_basic_set_get_dim(coef)));
909
910         coef = isl_basic_set_transform_dims(coef, isl_dim_set,
911                     isl_dim_size(dim, isl_dim_set), isl_mat_copy(src->cmap));
912         coef = isl_basic_set_transform_dims(coef, isl_dim_set,
913                     isl_dim_size(dim, isl_dim_set) + src->nvar,
914                     isl_mat_copy(dst->cmap));
915
916         nparam = isl_dim_size(src->dim, isl_dim_param);
917         total = isl_basic_set_total_dim(graph->lp);
918         dim_map = isl_dim_map_alloc(ctx, total);
919
920         isl_dim_map_range(dim_map, 1, 0, 0, 0, 1, 1);
921         isl_dim_map_range(dim_map, 4, 2, 1, 1, nparam, -1);
922         isl_dim_map_range(dim_map, 5, 2, 1, 1, nparam, 1);
923
924         isl_dim_map_range(dim_map, dst->start, 0, 0, 0, 1, -s);
925         isl_dim_map_range(dim_map, dst->start + 1, 2, 1, 1, dst->nparam, s);
926         isl_dim_map_range(dim_map, dst->start + 2, 2, 1, 1, dst->nparam, -s);
927         isl_dim_map_range(dim_map, dst->start + 2 * dst->nparam + 1, 2,
928                           isl_dim_size(dim, isl_dim_set) + src->nvar, 1,
929                           dst->nvar, s);
930         isl_dim_map_range(dim_map, dst->start + 2 * dst->nparam + 2, 2,
931                           isl_dim_size(dim, isl_dim_set) + src->nvar, 1,
932                           dst->nvar, -s);
933
934         isl_dim_map_range(dim_map, src->start, 0, 0, 0, 1, s);
935         isl_dim_map_range(dim_map, src->start + 1, 2, 1, 1, src->nparam, -s);
936         isl_dim_map_range(dim_map, src->start + 2, 2, 1, 1, src->nparam, s);
937         isl_dim_map_range(dim_map, src->start + 2 * src->nparam + 1, 2,
938                           isl_dim_size(dim, isl_dim_set), 1,
939                           src->nvar, -s);
940         isl_dim_map_range(dim_map, src->start + 2 * src->nparam + 2, 2,
941                           isl_dim_size(dim, isl_dim_set), 1,
942                           src->nvar, s);
943
944         graph->lp = isl_basic_set_extend_constraints(graph->lp,
945                         coef->n_eq, coef->n_ineq);
946         graph->lp = isl_basic_set_add_constraints_dim_map(graph->lp,
947                                                            coef, dim_map);
948         isl_dim_free(dim);
949
950         return 0;
951 }
952
953 static int add_all_validity_constraints(struct isl_sched_graph *graph)
954 {
955         int i;
956
957         for (i = 0; i < graph->n_edge; ++i) {
958                 struct isl_sched_edge *edge= &graph->edge[i];
959                 if (!edge->validity)
960                         continue;
961                 if (edge->src != edge->dst)
962                         continue;
963                 if (add_intra_validity_constraints(graph, edge) < 0)
964                         return -1;
965         }
966
967         for (i = 0; i < graph->n_edge; ++i) {
968                 struct isl_sched_edge *edge = &graph->edge[i];
969                 if (!edge->validity)
970                         continue;
971                 if (edge->src == edge->dst)
972                         continue;
973                 if (add_inter_validity_constraints(graph, edge) < 0)
974                         return -1;
975         }
976
977         return 0;
978 }
979
980 /* Add constraints to graph->lp that bound the dependence distance
981  * for all dependence relations.
982  * If a given proximity dependence is identical to a validity
983  * dependence, then the dependence distance is already bounded
984  * from below (by zero), so we only need to bound the distance
985  * from above.
986  * Otherwise, we need to bound the distance both from above and from below.
987  */
988 static int add_all_proximity_constraints(struct isl_sched_graph *graph)
989 {
990         int i;
991
992         for (i = 0; i < graph->n_edge; ++i) {
993                 struct isl_sched_edge *edge= &graph->edge[i];
994                 if (!edge->proximity)
995                         continue;
996                 if (edge->src == edge->dst &&
997                     add_intra_proximity_constraints(graph, edge, 1) < 0)
998                         return -1;
999                 if (edge->src != edge->dst &&
1000                     add_inter_proximity_constraints(graph, edge, 1) < 0)
1001                         return -1;
1002                 if (edge->validity)
1003                         continue;
1004                 if (edge->src == edge->dst &&
1005                     add_intra_proximity_constraints(graph, edge, -1) < 0)
1006                         return -1;
1007                 if (edge->src != edge->dst &&
1008                     add_inter_proximity_constraints(graph, edge, -1) < 0)
1009                         return -1;
1010         }
1011
1012         return 0;
1013 }
1014
1015 /* Compute a basis for the rows in the linear part of the schedule
1016  * and extend this basis to a full basis.  The remaining rows
1017  * can then be used to force linear independence from the rows
1018  * in the schedule.
1019  *
1020  * In particular, given the schedule rows S, we compute
1021  *
1022  *      S = H Q
1023  *
1024  * with H the Hermite normal form of S.  That is, all but the
1025  * first rank columns of Q are zero and so each row in S is
1026  * a linear combination of the first rank rows of Q.
1027  * The matrix Q is then transposed because we will write the
1028  * coefficients of the next schedule row as a column vector s
1029  * and express this s as a linear combination s = Q c of the
1030  * computed basis.
1031  */
1032 static int node_update_cmap(struct isl_sched_node *node)
1033 {
1034         isl_mat *H, *Q;
1035         int n_row = isl_mat_rows(node->sched);
1036
1037         H = isl_mat_sub_alloc(node->sched, 0, n_row,
1038                               1 + node->nparam, node->nvar);
1039
1040         H = isl_mat_left_hermite(H, 0, NULL, &Q);
1041         isl_mat_free(node->cmap);
1042         node->cmap = isl_mat_transpose(Q);
1043         node->rank = isl_mat_initial_non_zero_cols(H);
1044         isl_mat_free(H);
1045
1046         if (!node->cmap || node->rank < 0)
1047                 return -1;
1048         return 0;
1049 }
1050
1051 /* Count the number of equality and inequality constraints
1052  * that will be added.  If once is set, then we count
1053  * each edge exactly once.  Otherwise, we count as follows
1054  * validity             -> 1 (>= 0)
1055  * validity+proximity   -> 2 (>= 0 and upper bound)
1056  * proximity            -> 2 (lower and upper bound)
1057  */
1058 static int count_constraints(struct isl_sched_graph *graph,
1059         int *n_eq, int *n_ineq, int once)
1060 {
1061         int i;
1062         isl_basic_set *coef;
1063
1064         *n_eq = *n_ineq = 0;
1065         for (i = 0; i < graph->n_edge; ++i) {
1066                 struct isl_sched_edge *edge= &graph->edge[i];
1067                 isl_map *map = isl_map_copy(edge->map);
1068                 int f = once ? 1 : edge->proximity ? 2 : 1;
1069
1070                 if (edge->src == edge->dst)
1071                         coef = intra_coefficients(graph, map);
1072                 else
1073                         coef = inter_coefficients(graph, map);
1074                 if (!coef)
1075                         return -1;
1076                 *n_eq += f * coef->n_eq;
1077                 *n_ineq += f * coef->n_ineq;
1078                 isl_basic_set_free(coef);
1079         }
1080
1081         return 0;
1082 }
1083
1084 /* Construct an ILP problem for finding schedule coefficients
1085  * that result in non-negative, but small dependence distances
1086  * over all dependences.
1087  * In particular, the dependence distances over proximity edges
1088  * are bounded by m_0 + m_n n and we compute schedule coefficients
1089  * with small values (preferably zero) of m_n and m_0.
1090  *
1091  * All variables of the ILP are non-negative.  The actual coefficients
1092  * may be negative, so each coefficient is represented as the difference
1093  * of two non-negative variables.  The negative part always appears
1094  * immediately before the positive part.
1095  * Other than that, the variables have the following order
1096  *
1097  *      - sum of positive and negative parts of m_n coefficients
1098  *      - m_0
1099  *      - sum of positive and negative parts of all c_n coefficients
1100  *              (unconstrained when computing non-parametric schedules)
1101  *      - sum of positive and negative parts of all c_x coefficients
1102  *      - positive and negative parts of m_n coefficients
1103  *      - for each node
1104  *              - c_i_0
1105  *              - positive and negative parts of c_i_n (if parametric)
1106  *              - positive and negative parts of c_i_x
1107  *
1108  * The c_i_x are not represented directly, but through the columns of
1109  * node->cmap.  That is, the computed values are for variable t_i_x
1110  * such that c_i_x = Q t_i_x with Q equal to node->cmap.
1111  *
1112  * The constraints are those from the edges plus two or three equalities
1113  * to express the sums.
1114  */
1115 static int setup_lp(isl_ctx *ctx, struct isl_sched_graph *graph)
1116 {
1117         int i, j;
1118         int k;
1119         unsigned nparam;
1120         unsigned total;
1121         isl_dim *dim;
1122         int parametric;
1123         int param_pos;
1124         int n_eq, n_ineq;
1125
1126         parametric = ctx->opt->schedule_parametric;
1127         nparam = isl_dim_size(graph->node[0].dim, isl_dim_param);
1128         param_pos = 4;
1129         total = param_pos + 2 * nparam;
1130         for (i = 0; i < graph->n; ++i) {
1131                 struct isl_sched_node *node = &graph->node[graph->sorted[i]];
1132                 if (node_update_cmap(node) < 0)
1133                         return -1;
1134                 node->start = total;
1135                 total += 1 + 2 * (node->nparam + node->nvar);
1136         }
1137
1138         if (count_constraints(graph, &n_eq, &n_ineq, 0) < 0)
1139                 return -1;
1140
1141         dim = isl_dim_set_alloc(ctx, 0, total);
1142         isl_basic_set_free(graph->lp);
1143         n_eq += 2 + parametric;
1144         graph->lp = isl_basic_set_alloc_dim(dim, 0, n_eq, n_ineq);
1145
1146         k = isl_basic_set_alloc_equality(graph->lp);
1147         if (k < 0)
1148                 return -1;
1149         isl_seq_clr(graph->lp->eq[k], 1 +  total);
1150         isl_int_set_si(graph->lp->eq[k][1], -1);
1151         for (i = 0; i < 2 * nparam; ++i)
1152                 isl_int_set_si(graph->lp->eq[k][1 + param_pos + i], 1);
1153
1154         if (parametric) {
1155                 k = isl_basic_set_alloc_equality(graph->lp);
1156                 if (k < 0)
1157                         return -1;
1158                 isl_seq_clr(graph->lp->eq[k], 1 +  total);
1159                 isl_int_set_si(graph->lp->eq[k][3], -1);
1160                 for (i = 0; i < graph->n; ++i) {
1161                         int pos = 1 + graph->node[i].start + 1;
1162
1163                         for (j = 0; j < 2 * graph->node[i].nparam; ++j)
1164                                 isl_int_set_si(graph->lp->eq[k][pos + j], 1);
1165                 }
1166         }
1167
1168         k = isl_basic_set_alloc_equality(graph->lp);
1169         if (k < 0)
1170                 return -1;
1171         isl_seq_clr(graph->lp->eq[k], 1 +  total);
1172         isl_int_set_si(graph->lp->eq[k][4], -1);
1173         for (i = 0; i < graph->n; ++i) {
1174                 struct isl_sched_node *node = &graph->node[i];
1175                 int pos = 1 + node->start + 1 + 2 * node->nparam;
1176
1177                 for (j = 0; j < 2 * node->nvar; ++j)
1178                         isl_int_set_si(graph->lp->eq[k][pos + j], 1);
1179         }
1180
1181         if (add_all_validity_constraints(graph) < 0)
1182                 return -1;
1183         if (add_all_proximity_constraints(graph) < 0)
1184                 return -1;
1185
1186         return 0;
1187 }
1188
1189 /* Analyze the conflicting constraint found by
1190  * isl_tab_basic_set_non_trivial_lexmin.  If it corresponds to the validity
1191  * constraint of one of the edges between distinct nodes, living, moreover
1192  * in distinct SCCs, then record the source and sink SCC as this may
1193  * be a good place to cut between SCCs.
1194  */
1195 static int check_conflict(int con, void *user)
1196 {
1197         int i;
1198         struct isl_sched_graph *graph = user;
1199
1200         if (graph->src_scc >= 0)
1201                 return 0;
1202
1203         con -= graph->lp->n_eq;
1204
1205         if (con >= graph->lp->n_ineq)
1206                 return 0;
1207
1208         for (i = 0; i < graph->n_edge; ++i) {
1209                 if (!graph->edge[i].validity)
1210                         continue;
1211                 if (graph->edge[i].src == graph->edge[i].dst)
1212                         continue;
1213                 if (graph->edge[i].src->scc == graph->edge[i].dst->scc)
1214                         continue;
1215                 if (graph->edge[i].start > con)
1216                         continue;
1217                 if (graph->edge[i].end <= con)
1218                         continue;
1219                 graph->src_scc = graph->edge[i].src->scc;
1220                 graph->dst_scc = graph->edge[i].dst->scc;
1221         }
1222
1223         return 0;
1224 }
1225
1226 /* Check whether the next schedule row of the given node needs to be
1227  * non-trivial.  Lower-dimensional domains may have some trivial rows,
1228  * but as soon as the number of remaining required non-trivial rows
1229  * is as large as the number or remaining rows to be computed,
1230  * all remaining rows need to be non-trivial.
1231  */
1232 static int needs_row(struct isl_sched_graph *graph, struct isl_sched_node *node)
1233 {
1234         return node->nvar - node->rank >= graph->maxvar - graph->n_row;
1235 }
1236
1237 /* Solve the ILP problem constructed in setup_lp.
1238  * For each node such that all the remaining rows of its schedule
1239  * need to be non-trivial, we construct a non-triviality region.
1240  * This region imposes that the next row is independent of previous rows.
1241  * In particular the coefficients c_i_x are represented by t_i_x
1242  * variables with c_i_x = Q t_i_x and Q a unimodular matrix such that
1243  * its first columns span the rows of the previously computed part
1244  * of the schedule.  The non-triviality region enforces that at least
1245  * one of the remaining components of t_i_x is non-zero, i.e.,
1246  * that the new schedule row depends on at least one of the remaining
1247  * columns of Q.
1248  */
1249 static __isl_give isl_vec *solve_lp(struct isl_sched_graph *graph)
1250 {
1251         int i;
1252         isl_vec *sol;
1253         isl_basic_set *lp;
1254
1255         for (i = 0; i < graph->n; ++i) {
1256                 struct isl_sched_node *node = &graph->node[i];
1257                 int skip = node->rank;
1258                 graph->region[i].pos = node->start + 1 + 2*(node->nparam+skip);
1259                 if (needs_row(graph, node))
1260                         graph->region[i].len = 2 * (node->nvar - skip);
1261                 else
1262                         graph->region[i].len = 0;
1263         }
1264         lp = isl_basic_set_copy(graph->lp);
1265         sol = isl_tab_basic_set_non_trivial_lexmin(lp, 2, graph->n,
1266                                        graph->region, &check_conflict, graph);
1267         return sol;
1268 }
1269
1270 /* Update the schedules of all nodes based on the given solution
1271  * of the LP problem.
1272  * The new row is added to the current band.
1273  * All possibly negative coefficients are encoded as a difference
1274  * of two non-negative variables, so we need to perform the subtraction
1275  * here.  Moreover, if use_cmap is set, then the solution does
1276  * not refer to the actual coefficients c_i_x, but instead to variables
1277  * t_i_x such that c_i_x = Q t_i_x and Q is equal to node->cmap.
1278  * In this case, we then also need to perform this multiplication
1279  * to obtain the values of c_i_x.
1280  */
1281 static int update_schedule(struct isl_sched_graph *graph,
1282         __isl_take isl_vec *sol, int use_cmap)
1283 {
1284         int i, j;
1285         isl_vec *csol = NULL;
1286
1287         if (!sol)
1288                 goto error;
1289         if (sol->size == 0)
1290                 isl_die(sol->ctx, isl_error_internal,
1291                         "no solution found", goto error);
1292
1293         for (i = 0; i < graph->n; ++i) {
1294                 struct isl_sched_node *node = &graph->node[i];
1295                 int pos = node->start;
1296                 int row = isl_mat_rows(node->sched);
1297
1298                 isl_vec_free(csol);
1299                 csol = isl_vec_alloc(sol->ctx, node->nvar);
1300                 if (!csol)
1301                         goto error;
1302
1303                 isl_map_free(node->sched_map);
1304                 node->sched_map = NULL;
1305                 node->sched = isl_mat_add_rows(node->sched, 1);
1306                 if (!node->sched)
1307                         goto error;
1308                 node->sched = isl_mat_set_element(node->sched, row, 0,
1309                                                   sol->el[1 + pos]);
1310                 for (j = 0; j < node->nparam + node->nvar; ++j)
1311                         isl_int_sub(sol->el[1 + pos + 1 + 2 * j + 1],
1312                                     sol->el[1 + pos + 1 + 2 * j + 1],
1313                                     sol->el[1 + pos + 1 + 2 * j]);
1314                 for (j = 0; j < node->nparam; ++j)
1315                         node->sched = isl_mat_set_element(node->sched,
1316                                         row, 1 + j, sol->el[1+pos+1+2*j+1]);
1317                 for (j = 0; j < node->nvar; ++j)
1318                         isl_int_set(csol->el[j],
1319                                     sol->el[1+pos+1+2*(node->nparam+j)+1]);
1320                 if (use_cmap)
1321                         csol = isl_mat_vec_product(isl_mat_copy(node->cmap),
1322                                                    csol);
1323                 if (!csol)
1324                         goto error;
1325                 for (j = 0; j < node->nvar; ++j)
1326                         node->sched = isl_mat_set_element(node->sched,
1327                                         row, 1 + node->nparam + j, csol->el[j]);
1328                 node->band[graph->n_total_row] = graph->n_band;
1329         }
1330         isl_vec_free(sol);
1331         isl_vec_free(csol);
1332
1333         graph->n_row++;
1334         graph->n_total_row++;
1335
1336         return 0;
1337 error:
1338         isl_vec_free(sol);
1339         isl_vec_free(csol);
1340         return -1;
1341 }
1342
1343 /* Convert node->sched into a map and return this map.
1344  * We simply add equality constraints that express each output variable
1345  * as the affine combination of parameters and input variables specified
1346  * by the schedule matrix.
1347  *
1348  * The result is cached in node->sched_map, which needs to be released
1349  * whenever node->sched is updated.
1350  */
1351 static __isl_give isl_map *node_extract_schedule(struct isl_sched_node *node)
1352 {
1353         int i, j;
1354         isl_dim *dim;
1355         isl_basic_map *bmap;
1356         isl_constraint *c;
1357         int nrow, ncol;
1358         isl_int v;
1359
1360         if (node->sched_map)
1361                 return isl_map_copy(node->sched_map);
1362
1363         nrow = isl_mat_rows(node->sched);
1364         ncol = isl_mat_cols(node->sched) - 1;
1365         dim = isl_dim_from_domain(isl_dim_copy(node->dim));
1366         dim = isl_dim_add(dim, isl_dim_out, nrow);
1367         bmap = isl_basic_map_universe(isl_dim_copy(dim));
1368
1369         isl_int_init(v);
1370
1371         for (i = 0; i < nrow; ++i) {
1372                 c = isl_equality_alloc(isl_dim_copy(dim));
1373                 isl_constraint_set_coefficient_si(c, isl_dim_out, i, -1);
1374                 isl_mat_get_element(node->sched, i, 0, &v);
1375                 isl_constraint_set_constant(c, v);
1376                 for (j = 0; j < node->nparam; ++j) {
1377                         isl_mat_get_element(node->sched, i, 1 + j, &v);
1378                         isl_constraint_set_coefficient(c, isl_dim_param, j, v);
1379                 }
1380                 for (j = 0; j < node->nvar; ++j) {
1381                         isl_mat_get_element(node->sched,
1382                                             i, 1 + node->nparam + j, &v);
1383                         isl_constraint_set_coefficient(c, isl_dim_in, j, v);
1384                 }
1385                 bmap = isl_basic_map_add_constraint(bmap, c);
1386         }
1387
1388         isl_int_clear(v);
1389
1390         isl_dim_free(dim);
1391
1392         node->sched_map = isl_map_from_basic_map(bmap);
1393         return isl_map_copy(node->sched_map);
1394 }
1395
1396 /* Update the given dependence relation based on the current schedule.
1397  * That is, intersect the dependence relation with a map expressing
1398  * that source and sink are executed within the same iteration of
1399  * the current schedule.
1400  * This is not the most efficient way, but this shouldn't be a critical
1401  * operation.
1402  */
1403 static __isl_give isl_map *specialize(__isl_take isl_map *map,
1404         struct isl_sched_node *src, struct isl_sched_node *dst)
1405 {
1406         isl_map *src_sched, *dst_sched, *id;
1407
1408         src_sched = node_extract_schedule(src);
1409         dst_sched = node_extract_schedule(dst);
1410         id = isl_map_apply_range(src_sched, isl_map_reverse(dst_sched));
1411         return isl_map_intersect(map, id);
1412 }
1413
1414 /* Update the dependence relations of all edges based on the current schedule.
1415  * If a dependence is carried completely by the current schedule, then
1416  * it is removed and edge_table is updated accordingly.
1417  */
1418 static int update_edges(isl_ctx *ctx, struct isl_sched_graph *graph)
1419 {
1420         int i;
1421         int reset_table = 0;
1422
1423         for (i = graph->n_edge - 1; i >= 0; --i) {
1424                 struct isl_sched_edge *edge = &graph->edge[i];
1425                 edge->map = specialize(edge->map, edge->src, edge->dst);
1426                 if (!edge->map)
1427                         return -1;
1428
1429                 if (isl_map_fast_is_empty(edge->map)) {
1430                         reset_table = 1;
1431                         isl_map_free(edge->map);
1432                         if (i != graph->n_edge - 1)
1433                                 graph->edge[i] = graph->edge[graph->n_edge - 1];
1434                         graph->n_edge--;
1435                 }
1436         }
1437
1438         if (reset_table) {
1439                 isl_hash_table_free(ctx, graph->edge_table);
1440                 graph->edge_table = NULL;
1441                 return graph_init_edge_table(ctx, graph);
1442         }
1443
1444         return 0;
1445 }
1446
1447 static void next_band(struct isl_sched_graph *graph)
1448 {
1449         graph->band_start = graph->n_total_row;
1450         graph->n_band++;
1451 }
1452
1453 /* Topologically sort statements mapped to same schedule iteration
1454  * and add a row to the schedule corresponding to this order.
1455  */
1456 static int sort_statements(isl_ctx *ctx, struct isl_sched_graph *graph)
1457 {
1458         int i, j;
1459
1460         if (graph->n <= 1)
1461                 return 0;
1462
1463         if (update_edges(ctx, graph) < 0)
1464                 return -1;
1465
1466         if (graph->n_edge == 0)
1467                 return 0;
1468
1469         if (detect_sccs(graph) < 0)
1470                 return -1;
1471
1472         for (i = 0; i < graph->n; ++i) {
1473                 struct isl_sched_node *node = &graph->node[i];
1474                 int row = isl_mat_rows(node->sched);
1475                 int cols = isl_mat_cols(node->sched);
1476
1477                 isl_map_free(node->sched_map);
1478                 node->sched_map = NULL;
1479                 node->sched = isl_mat_add_rows(node->sched, 1);
1480                 if (!node->sched)
1481                         return -1;
1482                 node->sched = isl_mat_set_element_si(node->sched, row, 0,
1483                                                      node->scc);
1484                 for (j = 1; j < cols; ++j)
1485                         node->sched = isl_mat_set_element_si(node->sched,
1486                                                              row, j, 0);
1487                 node->band[graph->n_total_row] = graph->n_band;
1488         }
1489
1490         graph->n_total_row++;
1491         next_band(graph);
1492
1493         return 0;
1494 }
1495
1496 /* Construct an isl_schedule based on the computed schedule stored
1497  * in graph and with parameters specified by dim.
1498  */
1499 static __isl_give isl_schedule *extract_schedule(struct isl_sched_graph *graph,
1500         __isl_take isl_dim *dim)
1501 {
1502         int i;
1503         isl_ctx *ctx;
1504         isl_schedule *sched = NULL;
1505                 
1506         if (!dim)
1507                 return NULL;
1508
1509         ctx = isl_dim_get_ctx(dim);
1510         sched = isl_calloc(ctx, struct isl_schedule,
1511                            sizeof(struct isl_schedule) +
1512                            (graph->n - 1) * sizeof(struct isl_schedule_node));
1513         if (!sched)
1514                 goto error;
1515
1516         sched->n = graph->n;
1517         sched->n_band = graph->n_band;
1518         sched->n_total_row = graph->n_total_row;
1519
1520         for (i = 0; i < sched->n; ++i) {
1521                 int r, b;
1522                 int *band_end;
1523
1524                 band_end = isl_alloc_array(ctx, int, graph->n_band);
1525                 if (!band_end)
1526                         goto error;
1527                 sched->node[i].sched = node_extract_schedule(&graph->node[i]);
1528                 sched->node[i].band_end = band_end;
1529
1530                 for (r = b = 0; r < graph->n_total_row; ++r) {
1531                         if (graph->node[i].band[r] == b)
1532                                 continue;
1533                         band_end[b++] = r;
1534                         if (graph->node[i].band[r] == -1)
1535                                 break;
1536                 }
1537                 if (r == graph->n_total_row)
1538                         band_end[b++] = r;
1539                 sched->node[i].n_band = b;
1540         }
1541
1542         sched->dim = dim;
1543
1544         return sched;
1545 error:
1546         isl_dim_free(dim);
1547         isl_schedule_free(sched);
1548         return NULL;
1549 }
1550
1551 /* Copy nodes that satisfy node_pred from the src dependence graph
1552  * to the dst dependence graph.
1553  */
1554 static int copy_nodes(struct isl_sched_graph *dst, struct isl_sched_graph *src,
1555         int (*node_pred)(struct isl_sched_node *node, int data), int data)
1556 {
1557         int i;
1558
1559         dst->n = 0;
1560         for (i = 0; i < src->n; ++i) {
1561                 if (!node_pred(&src->node[i], data))
1562                         continue;
1563                 dst->node[dst->n].dim = isl_dim_copy(src->node[i].dim);
1564                 dst->node[dst->n].nvar = src->node[i].nvar;
1565                 dst->node[dst->n].nparam = src->node[i].nparam;
1566                 dst->node[dst->n].sched = isl_mat_copy(src->node[i].sched);
1567                 dst->node[dst->n].sched_map =
1568                         isl_map_copy(src->node[i].sched_map);
1569                 dst->node[dst->n].band = src->node[i].band;
1570                 dst->n++;
1571         }
1572
1573         return 0;
1574 }
1575
1576 /* Copy non-empty edges that satisfy edge_pred from the src dependence graph
1577  * to the dst dependence graph.
1578  */
1579 static int copy_edges(isl_ctx *ctx, struct isl_sched_graph *dst,
1580         struct isl_sched_graph *src,
1581         int (*edge_pred)(struct isl_sched_edge *edge, int data), int data)
1582 {
1583         int i;
1584
1585         dst->n_edge = 0;
1586         for (i = 0; i < src->n_edge; ++i) {
1587                 struct isl_sched_edge *edge = &src->edge[i];
1588                 isl_map *map;
1589
1590                 if (!edge_pred(edge, data))
1591                         continue;
1592
1593                 if (isl_map_fast_is_empty(edge->map))
1594                         continue;
1595
1596                 map = isl_map_copy(edge->map);
1597
1598                 dst->edge[dst->n_edge].src =
1599                         graph_find_node(ctx, dst, edge->src->dim);
1600                 dst->edge[dst->n_edge].dst =
1601                         graph_find_node(ctx, dst, edge->dst->dim);
1602                 dst->edge[dst->n_edge].map = map;
1603                 dst->edge[dst->n_edge].validity = edge->validity;
1604                 dst->edge[dst->n_edge].proximity = edge->proximity;
1605                 dst->n_edge++;
1606         }
1607
1608         return 0;
1609 }
1610
1611 /* Given a "src" dependence graph that contains the nodes from "dst"
1612  * that satisfy node_pred, copy the schedule computed in "src"
1613  * for those nodes back to "dst".
1614  */
1615 static int copy_schedule(struct isl_sched_graph *dst,
1616         struct isl_sched_graph *src,
1617         int (*node_pred)(struct isl_sched_node *node, int data), int data)
1618 {
1619         int i;
1620
1621         src->n = 0;
1622         for (i = 0; i < dst->n; ++i) {
1623                 if (!node_pred(&dst->node[i], data))
1624                         continue;
1625                 isl_mat_free(dst->node[i].sched);
1626                 isl_map_free(dst->node[i].sched_map);
1627                 dst->node[i].sched = isl_mat_copy(src->node[src->n].sched);
1628                 dst->node[i].sched_map =
1629                         isl_map_copy(src->node[src->n].sched_map);
1630                 src->n++;
1631         }
1632
1633         dst->n_total_row = src->n_total_row;
1634         dst->n_band = src->n_band;
1635
1636         return 0;
1637 }
1638
1639 /* Compute the maximal number of variables over all nodes.
1640  * This is the maximal number of linearly independent schedule
1641  * rows that we need to compute.
1642  * Just in case we end up in a part of the dependence graph
1643  * with only lower-dimensional domains, we make sure we will
1644  * compute the required amount of extra linearly independent rows.
1645  */
1646 static int compute_maxvar(struct isl_sched_graph *graph)
1647 {
1648         int i;
1649
1650         graph->maxvar = 0;
1651         for (i = 0; i < graph->n; ++i) {
1652                 struct isl_sched_node *node = &graph->node[i];
1653                 int nvar;
1654
1655                 if (node_update_cmap(node) < 0)
1656                         return -1;
1657                 nvar = node->nvar + graph->n_row - node->rank;
1658                 if (nvar > graph->maxvar)
1659                         graph->maxvar = nvar;
1660         }
1661
1662         return 0;
1663 }
1664
1665 static int compute_schedule(isl_ctx *ctx, struct isl_sched_graph *graph);
1666 static int compute_schedule_wcc(isl_ctx *ctx, struct isl_sched_graph *graph);
1667
1668 /* Compute a schedule for a subgraph of "graph".  In particular, for
1669  * the graph composed of nodes that satisfy node_pred and edges that
1670  * that satisfy edge_pred.  The caller should precompute the number
1671  * of nodes and edges that satisfy these predicates and pass them along
1672  * as "n" and "n_edge".
1673  * If the subgraph is known to consist of a single component, then wcc should
1674  * be set and then we call compute_schedule_wcc on the constructed subgraph.
1675  * Otherwise, we call compute_schedule, which will check whether the subgraph
1676  * is connected.
1677  */
1678 static int compute_sub_schedule(isl_ctx *ctx,
1679         struct isl_sched_graph *graph, int n, int n_edge,
1680         int (*node_pred)(struct isl_sched_node *node, int data),
1681         int (*edge_pred)(struct isl_sched_edge *edge, int data),
1682         int data, int wcc)
1683 {
1684         struct isl_sched_graph split = { 0 };
1685
1686         if (graph_alloc(ctx, &split, n, n_edge) < 0)
1687                 goto error;
1688         if (copy_nodes(&split, graph, node_pred, data) < 0)
1689                 goto error;
1690         if (graph_init_table(ctx, &split) < 0)
1691                 goto error;
1692         if (copy_edges(ctx, &split, graph, edge_pred, data) < 0)
1693                 goto error;
1694         if (graph_init_edge_table(ctx, &split) < 0)
1695                 goto error;
1696         split.n_row = graph->n_row;
1697         split.n_total_row = graph->n_total_row;
1698         split.n_band = graph->n_band;
1699         split.band_start = graph->band_start;
1700
1701         if (wcc && compute_schedule_wcc(ctx, &split) < 0)
1702                 goto error;
1703         if (!wcc && compute_schedule(ctx, &split) < 0)
1704                 goto error;
1705
1706         copy_schedule(graph, &split, node_pred, data);
1707
1708         graph_free(ctx, &split);
1709         return 0;
1710 error:
1711         graph_free(ctx, &split);
1712         return -1;
1713 }
1714
1715 static int node_scc_exactly(struct isl_sched_node *node, int scc)
1716 {
1717         return node->scc == scc;
1718 }
1719
1720 static int node_scc_at_most(struct isl_sched_node *node, int scc)
1721 {
1722         return node->scc <= scc;
1723 }
1724
1725 static int node_scc_at_least(struct isl_sched_node *node, int scc)
1726 {
1727         return node->scc >= scc;
1728 }
1729
1730 static int edge_src_scc_exactly(struct isl_sched_edge *edge, int scc)
1731 {
1732         return edge->src->scc == scc;
1733 }
1734
1735 static int edge_dst_scc_at_most(struct isl_sched_edge *edge, int scc)
1736 {
1737         return edge->dst->scc <= scc;
1738 }
1739
1740 static int edge_src_scc_at_least(struct isl_sched_edge *edge, int scc)
1741 {
1742         return edge->src->scc >= scc;
1743 }
1744
1745 /* Pad the schedules of all nodes with zero rows such that in the end
1746  * they all have graph->n_total_row rows.
1747  * The extra rows don't belong to any band, so they get assigned band number -1.
1748  */
1749 static int pad_schedule(struct isl_sched_graph *graph)
1750 {
1751         int i, j;
1752
1753         for (i = 0; i < graph->n; ++i) {
1754                 struct isl_sched_node *node = &graph->node[i];
1755                 int row = isl_mat_rows(node->sched);
1756                 if (graph->n_total_row > row) {
1757                         isl_map_free(node->sched_map);
1758                         node->sched_map = NULL;
1759                 }
1760                 node->sched = isl_mat_add_zero_rows(node->sched,
1761                                                     graph->n_total_row - row);
1762                 if (!node->sched)
1763                         return -1;
1764                 for (j = row; j < graph->n_total_row; ++j)
1765                         node->band[j] = -1;
1766         }
1767
1768         return 0;
1769 }
1770
1771 /* Split the current graph into two parts and compute a schedule for each
1772  * part individually.  In particular, one part consists of all SCCs up
1773  * to and including graph->src_scc, while the other part contains the other
1774  * SCCS.
1775  *
1776  * The split is enforced in the schedule by constant rows with two different
1777  * values (0 and 1).  These constant rows replace the previously computed rows
1778  * in the current band.
1779  * It would be possible to reuse them as the first rows in the next
1780  * band, but recomputing them may result in better rows as we are looking
1781  * at a smaller part of the dependence graph.
1782  */
1783 static int compute_split_schedule(isl_ctx *ctx, struct isl_sched_graph *graph)
1784 {
1785         int i, j, n, e1, e2;
1786         int n_total_row, orig_total_row;
1787         int n_band, orig_band;
1788         int drop;
1789
1790         drop = graph->n_total_row - graph->band_start;
1791         graph->n_total_row -= drop;
1792         graph->n_row -= drop;
1793
1794         n = 0;
1795         for (i = 0; i < graph->n; ++i) {
1796                 struct isl_sched_node *node = &graph->node[i];
1797                 int row = isl_mat_rows(node->sched) - drop;
1798                 int cols = isl_mat_cols(node->sched);
1799                 int before = node->scc <= graph->src_scc;
1800
1801                 if (before)
1802                         n++;
1803
1804                 isl_map_free(node->sched_map);
1805                 node->sched_map = NULL;
1806                 node->sched = isl_mat_drop_rows(node->sched,
1807                                                 graph->band_start, drop);
1808                 node->sched = isl_mat_add_rows(node->sched, 1);
1809                 if (!node->sched)
1810                         return -1;
1811                 node->sched = isl_mat_set_element_si(node->sched, row, 0,
1812                                                      !before);
1813                 for (j = 1; j < cols; ++j)
1814                         node->sched = isl_mat_set_element_si(node->sched,
1815                                                              row, j, 0);
1816                 node->band[graph->n_total_row] = graph->n_band;
1817         }
1818
1819         e1 = e2 = 0;
1820         for (i = 0; i < graph->n_edge; ++i) {
1821                 if (graph->edge[i].dst->scc <= graph->src_scc)
1822                         e1++;
1823                 if (graph->edge[i].src->scc > graph->src_scc)
1824                         e2++;
1825         }
1826
1827         graph->n_total_row++;
1828         next_band(graph);
1829
1830         orig_total_row = graph->n_total_row;
1831         orig_band = graph->n_band;
1832         if (compute_sub_schedule(ctx, graph, n, e1,
1833                                 &node_scc_at_most, &edge_dst_scc_at_most,
1834                                 graph->src_scc, 0) < 0)
1835                 return -1;
1836         n_total_row = graph->n_total_row;
1837         graph->n_total_row = orig_total_row;
1838         n_band = graph->n_band;
1839         graph->n_band = orig_band;
1840         if (compute_sub_schedule(ctx, graph, graph->n - n, e2,
1841                                 &node_scc_at_least, &edge_src_scc_at_least,
1842                                 graph->src_scc + 1, 0) < 0)
1843                 return -1;
1844         if (n_total_row > graph->n_total_row)
1845                 graph->n_total_row = n_total_row;
1846         if (n_band > graph->n_band)
1847                 graph->n_band = n_band;
1848
1849         return pad_schedule(graph);
1850 }
1851
1852 /* Compute the next band of the schedule after updating the dependence
1853  * relations based on the the current schedule.
1854  */
1855 static int compute_next_band(isl_ctx *ctx, struct isl_sched_graph *graph)
1856 {
1857         if (update_edges(ctx, graph) < 0)
1858                 return -1;
1859         next_band(graph);
1860                 
1861         return compute_schedule(ctx, graph);
1862 }
1863
1864 /* Add constraints to graph->lp that force the dependence of edge i
1865  * to be respected and attempt to carry it, where edge i is one from
1866  * a node j to itself.
1867  * That is, add constraints that enforce
1868  *
1869  *      (c_j_0 + c_j_n n + c_j_x y) - (c_j_0 + c_j_n n + c_j_x x)
1870  *      = c_j_x (y - x) >= e_i
1871  *
1872  * for each (x,y) in R.
1873  * We obtain general constraints on coefficients (c_0, c_n, c_x)
1874  * of valid constraints for (y - x) and then plug in (-e_i, 0, c_j_x),
1875  * with each coefficient in c_j_x represented as a pair of non-negative
1876  * coefficients.
1877  */
1878 static int add_intra_constraints(struct isl_sched_graph *graph, int i)
1879 {
1880         unsigned total;
1881         struct isl_sched_edge *edge= &graph->edge[i];
1882         isl_map *map = isl_map_copy(edge->map);
1883         isl_ctx *ctx = isl_map_get_ctx(map);
1884         isl_dim *dim;
1885         isl_dim_map *dim_map;
1886         isl_basic_set *coef;
1887         struct isl_sched_node *node = edge->src;
1888
1889         coef = intra_coefficients(graph, map);
1890
1891         dim = isl_dim_domain(isl_dim_unwrap(isl_basic_set_get_dim(coef)));
1892
1893         total = isl_basic_set_total_dim(graph->lp);
1894         dim_map = isl_dim_map_alloc(ctx, total);
1895         isl_dim_map_range(dim_map, 3 + i, 0, 0, 0, 1, -1);
1896         isl_dim_map_range(dim_map, node->start + 2 * node->nparam + 1, 2,
1897                           isl_dim_size(dim, isl_dim_set), 1,
1898                           node->nvar, -1);
1899         isl_dim_map_range(dim_map, node->start + 2 * node->nparam + 2, 2,
1900                           isl_dim_size(dim, isl_dim_set), 1,
1901                           node->nvar, 1);
1902         graph->lp = isl_basic_set_extend_constraints(graph->lp,
1903                         coef->n_eq, coef->n_ineq);
1904         graph->lp = isl_basic_set_add_constraints_dim_map(graph->lp,
1905                                                            coef, dim_map);
1906         isl_dim_free(dim);
1907
1908         return 0;
1909 }
1910
1911 /* Add constraints to graph->lp that force the dependence of edge i
1912  * to be respected and attempt to carry it, where edge i is one from
1913  * node j to node k.
1914  * That is, add constraints that enforce
1915  *
1916  *      (c_k_0 + c_k_n n + c_k_x y) - (c_j_0 + c_j_n n + c_j_x x) >= e_i
1917  *
1918  * for each (x,y) in R.
1919  * We obtain general constraints on coefficients (c_0, c_n, c_x)
1920  * of valid constraints for R and then plug in
1921  * (-e_i + c_k_0 - c_j_0, c_k_n - c_j_n, c_k_x - c_j_x)
1922  * with each coefficient (except e_i, c_k_0 and c_j_0)
1923  * represented as a pair of non-negative coefficients.
1924  */
1925 static int add_inter_constraints(struct isl_sched_graph *graph, int i)
1926 {
1927         unsigned total;
1928         struct isl_sched_edge *edge= &graph->edge[i];
1929         isl_map *map = isl_map_copy(edge->map);
1930         isl_ctx *ctx = isl_map_get_ctx(map);
1931         isl_dim *dim;
1932         isl_dim_map *dim_map;
1933         isl_basic_set *coef;
1934         struct isl_sched_node *src = edge->src;
1935         struct isl_sched_node *dst = edge->dst;
1936
1937         coef = inter_coefficients(graph, map);
1938
1939         dim = isl_dim_domain(isl_dim_unwrap(isl_basic_set_get_dim(coef)));
1940
1941         total = isl_basic_set_total_dim(graph->lp);
1942         dim_map = isl_dim_map_alloc(ctx, total);
1943
1944         isl_dim_map_range(dim_map, 3 + i, 0, 0, 0, 1, -1);
1945
1946         isl_dim_map_range(dim_map, dst->start, 0, 0, 0, 1, 1);
1947         isl_dim_map_range(dim_map, dst->start + 1, 2, 1, 1, dst->nparam, -1);
1948         isl_dim_map_range(dim_map, dst->start + 2, 2, 1, 1, dst->nparam, 1);
1949         isl_dim_map_range(dim_map, dst->start + 2 * dst->nparam + 1, 2,
1950                           isl_dim_size(dim, isl_dim_set) + src->nvar, 1,
1951                           dst->nvar, -1);
1952         isl_dim_map_range(dim_map, dst->start + 2 * dst->nparam + 2, 2,
1953                           isl_dim_size(dim, isl_dim_set) + src->nvar, 1,
1954                           dst->nvar, 1);
1955
1956         isl_dim_map_range(dim_map, src->start, 0, 0, 0, 1, -1);
1957         isl_dim_map_range(dim_map, src->start + 1, 2, 1, 1, src->nparam, 1);
1958         isl_dim_map_range(dim_map, src->start + 2, 2, 1, 1, src->nparam, -1);
1959         isl_dim_map_range(dim_map, src->start + 2 * src->nparam + 1, 2,
1960                           isl_dim_size(dim, isl_dim_set), 1,
1961                           src->nvar, 1);
1962         isl_dim_map_range(dim_map, src->start + 2 * src->nparam + 2, 2,
1963                           isl_dim_size(dim, isl_dim_set), 1,
1964                           src->nvar, -1);
1965
1966         graph->lp = isl_basic_set_extend_constraints(graph->lp,
1967                         coef->n_eq, coef->n_ineq);
1968         graph->lp = isl_basic_set_add_constraints_dim_map(graph->lp,
1969                                                            coef, dim_map);
1970         isl_dim_free(dim);
1971
1972         return 0;
1973 }
1974
1975 /* Add constraints to graph->lp that force all dependence
1976  * to be respected and attempt to carry it.
1977  */
1978 static int add_all_constraints(struct isl_sched_graph *graph)
1979 {
1980         int i;
1981
1982         for (i = 0; i < graph->n_edge; ++i) {
1983                 struct isl_sched_edge *edge= &graph->edge[i];
1984                 if (edge->src == edge->dst &&
1985                     add_intra_constraints(graph, i) < 0)
1986                         return -1;
1987                 if (edge->src != edge->dst &&
1988                     add_inter_constraints(graph, i) < 0)
1989                         return -1;
1990         }
1991
1992         return 0;
1993 }
1994
1995 /* Construct an LP problem for finding schedule coefficients
1996  * such that the schedule carries as many dependences as possible.
1997  * In particular, for each dependence i, we bound the dependence distance
1998  * from below by e_i, with 0 <= e_i <= 1 and then maximize the sum
1999  * of all e_i's.  Dependence with e_i = 0 in the solution are simply
2000  * respected, while those with e_i > 0 (in practice e_i = 1) are carried.
2001  *
2002  * All variables of the LP are non-negative.  The actual coefficients
2003  * may be negative, so each coefficient is represented as the difference
2004  * of two non-negative variables.  The negative part always appears
2005  * immediately before the positive part.
2006  * Other than that, the variables have the following order
2007  *
2008  *      - sum of (1 - e_i) over all edges
2009  *      - sum of positive and negative parts of all c_n coefficients
2010  *              (unconstrained when computing non-parametric schedules)
2011  *      - sum of positive and negative parts of all c_x coefficients
2012  *      - for each edge
2013  *              - e_i
2014  *      - for each node
2015  *              - c_i_0
2016  *              - positive and negative parts of c_i_n (if parametric)
2017  *              - positive and negative parts of c_i_x
2018  *
2019  * The constraints are those from the edges plus three equalities
2020  * to express the sums and n_edge inequalities to express e_i <= 1.
2021  */
2022 static int setup_carry_lp(isl_ctx *ctx, struct isl_sched_graph *graph)
2023 {
2024         int i, j;
2025         int k;
2026         isl_dim *dim;
2027         unsigned total;
2028         int n_eq, n_ineq;
2029
2030         total = 3 + graph->n_edge;
2031         for (i = 0; i < graph->n; ++i) {
2032                 struct isl_sched_node *node = &graph->node[graph->sorted[i]];
2033                 node->start = total;
2034                 total += 1 + 2 * (node->nparam + node->nvar);
2035         }
2036
2037         if (count_constraints(graph, &n_eq, &n_ineq, 1) < 0)
2038                 return -1;
2039
2040         dim = isl_dim_set_alloc(ctx, 0, total);
2041         isl_basic_set_free(graph->lp);
2042         n_eq += 3;
2043         n_ineq += graph->n_edge;
2044         graph->lp = isl_basic_set_alloc_dim(dim, 0, n_eq, n_ineq);
2045         graph->lp = isl_basic_set_set_rational(graph->lp);
2046
2047         k = isl_basic_set_alloc_equality(graph->lp);
2048         if (k < 0)
2049                 return -1;
2050         isl_seq_clr(graph->lp->eq[k], 1 +  total);
2051         isl_int_set_si(graph->lp->eq[k][0], -graph->n_edge);
2052         isl_int_set_si(graph->lp->eq[k][1], 1);
2053         for (i = 0; i < graph->n_edge; ++i)
2054                 isl_int_set_si(graph->lp->eq[k][4 + i], 1);
2055
2056         k = isl_basic_set_alloc_equality(graph->lp);
2057         if (k < 0)
2058                 return -1;
2059         isl_seq_clr(graph->lp->eq[k], 1 +  total);
2060         isl_int_set_si(graph->lp->eq[k][2], -1);
2061         for (i = 0; i < graph->n; ++i) {
2062                 int pos = 1 + graph->node[i].start + 1;
2063
2064                 for (j = 0; j < 2 * graph->node[i].nparam; ++j)
2065                         isl_int_set_si(graph->lp->eq[k][pos + j], 1);
2066         }
2067
2068         k = isl_basic_set_alloc_equality(graph->lp);
2069         if (k < 0)
2070                 return -1;
2071         isl_seq_clr(graph->lp->eq[k], 1 +  total);
2072         isl_int_set_si(graph->lp->eq[k][3], -1);
2073         for (i = 0; i < graph->n; ++i) {
2074                 struct isl_sched_node *node = &graph->node[i];
2075                 int pos = 1 + node->start + 1 + 2 * node->nparam;
2076
2077                 for (j = 0; j < 2 * node->nvar; ++j)
2078                         isl_int_set_si(graph->lp->eq[k][pos + j], 1);
2079         }
2080
2081         for (i = 0; i < graph->n_edge; ++i) {
2082                 k = isl_basic_set_alloc_inequality(graph->lp);
2083                 if (k < 0)
2084                         return -1;
2085                 isl_seq_clr(graph->lp->ineq[k], 1 +  total);
2086                 isl_int_set_si(graph->lp->ineq[k][4 + i], -1);
2087                 isl_int_set_si(graph->lp->ineq[k][0], 1);
2088         }
2089
2090         if (add_all_constraints(graph) < 0)
2091                 return -1;
2092
2093         return 0;
2094 }
2095
2096 /* Construct a schedule row for each node such that as many dependences
2097  * as possible are carried and then continue with the next band.
2098  */
2099 static int carry_dependences(isl_ctx *ctx, struct isl_sched_graph *graph)
2100 {
2101         isl_vec *sol;
2102         isl_basic_set *lp;
2103
2104         if (setup_carry_lp(ctx, graph) < 0)
2105                 return -1;
2106
2107         lp = isl_basic_set_copy(graph->lp);
2108         sol = isl_tab_basic_set_non_neg_lexmin(lp);
2109         if (!sol)
2110                 return -1;
2111
2112         if (sol->size == 0) {
2113                 isl_vec_free(sol);
2114                 isl_die(ctx, isl_error_internal,
2115                         "error in schedule construction", return -1);
2116         }
2117
2118         if (isl_int_cmp_si(sol->el[1], graph->n_edge) >= 0) {
2119                 isl_vec_free(sol);
2120                 isl_die(ctx, isl_error_unknown,
2121                         "unable to carry dependences", return -1);
2122         }
2123
2124         if (update_schedule(graph, sol, 0) < 0)
2125                 return -1;
2126
2127         return compute_next_band(ctx, graph);
2128 }
2129
2130 /* Compute a schedule for a connected dependence graph.
2131  * We try to find a sequence of as many schedule rows as possible that result
2132  * in non-negative dependence distances (independent of the previous rows
2133  * in the sequence, i.e., such that the sequence is tilable).
2134  * If we can't find any more rows we either
2135  * - split between SCCs and start over (assuming we found an interesting
2136  *      pair of SCCs between which to split)
2137  * - continue with the next band (assuming the current band has at least
2138  *      one row)
2139  * - try to carry as many dependences as possible and continue with the next
2140  *      band
2141  *
2142  * If we manage to complete the schedule, we finish off by topologically
2143  * sorting the statements based on the remaining dependences.
2144  */
2145 static int compute_schedule_wcc(isl_ctx *ctx, struct isl_sched_graph *graph)
2146 {
2147         if (detect_sccs(graph) < 0)
2148                 return -1;
2149         sort_sccs(graph);
2150
2151         if (compute_maxvar(graph) < 0)
2152                 return -1;
2153
2154         while (graph->n_row < graph->maxvar) {
2155                 isl_vec *sol;
2156
2157                 graph->src_scc = -1;
2158                 graph->dst_scc = -1;
2159
2160                 if (setup_lp(ctx, graph) < 0)
2161                         return -1;
2162                 sol = solve_lp(graph);
2163                 if (!sol)
2164                         return -1;
2165                 if (sol->size == 0) {
2166                         isl_vec_free(sol);
2167                         if (graph->src_scc >= 0)
2168                                 return compute_split_schedule(ctx, graph);
2169                         if (graph->n_total_row > graph->band_start)
2170                                 return compute_next_band(ctx, graph);
2171                         return carry_dependences(ctx, graph);
2172                 }
2173                 if (update_schedule(graph, sol, 1) < 0)
2174                         return -1;
2175         }
2176
2177         if (graph->n_total_row > graph->band_start)
2178                 next_band(graph);
2179         return sort_statements(ctx, graph);
2180 }
2181
2182 /* Compute a schedule for each component (identified by node->scc)
2183  * of the dependence graph separately and then combine the results.
2184  */
2185 static int compute_component_schedule(isl_ctx *ctx,
2186         struct isl_sched_graph *graph)
2187 {
2188         int wcc, i;
2189         int n, n_edge;
2190         int n_total_row, orig_total_row;
2191         int n_band, orig_band;
2192
2193         n_total_row = 0;
2194         orig_total_row = graph->n_total_row;
2195         n_band = 0;
2196         orig_band = graph->n_band;
2197         for (wcc = 0; wcc < graph->scc; ++wcc) {
2198                 n = 0;
2199                 for (i = 0; i < graph->n; ++i)
2200                         if (graph->node[i].scc == wcc)
2201                                 n++;
2202                 n_edge = 0;
2203                 for (i = 0; i < graph->n_edge; ++i)
2204                         if (graph->edge[i].src->scc == wcc)
2205                                 n_edge++;
2206
2207                 if (compute_sub_schedule(ctx, graph, n, n_edge,
2208                                     &node_scc_exactly,
2209                                     &edge_src_scc_exactly, wcc, 1) < 0)
2210                         return -1;
2211                 if (graph->n_total_row > n_total_row)
2212                         n_total_row = graph->n_total_row;
2213                 graph->n_total_row = orig_total_row;
2214                 if (graph->n_band > n_band)
2215                         n_band = graph->n_band;
2216                 graph->n_band = orig_band;
2217         }
2218
2219         graph->n_total_row = n_total_row;
2220         graph->n_band = n_band;
2221
2222         return pad_schedule(graph);
2223 }
2224
2225 /* Compute a schedule for the given dependence graph.
2226  * We first check if the graph is connected (through validity dependences)
2227  * and if so compute a schedule for each component separately.
2228  */
2229 static int compute_schedule(isl_ctx *ctx, struct isl_sched_graph *graph)
2230 {
2231         if (detect_wccs(graph) < 0)
2232                 return -1;
2233
2234         if (graph->scc > 1)
2235                 return compute_component_schedule(ctx, graph);
2236
2237         return compute_schedule_wcc(ctx, graph);
2238 }
2239
2240 /* Compute a schedule for the given union of domains that respects
2241  * all the validity dependences and tries to minimize the dependence
2242  * distances over the proximity dependences.
2243  */
2244 __isl_give isl_schedule *isl_union_set_compute_schedule(
2245         __isl_take isl_union_set *domain,
2246         __isl_take isl_union_map *validity,
2247         __isl_take isl_union_map *proximity)
2248 {
2249         isl_ctx *ctx = isl_union_set_get_ctx(domain);
2250         isl_dim *dim;
2251         struct isl_sched_graph graph = { 0 };
2252         isl_schedule *sched;
2253
2254         domain = isl_union_set_align_params(domain,
2255                                             isl_union_map_get_dim(validity));
2256         domain = isl_union_set_align_params(domain,
2257                                             isl_union_map_get_dim(proximity));
2258         dim = isl_union_set_get_dim(domain);
2259         validity = isl_union_map_align_params(validity, isl_dim_copy(dim));
2260         proximity = isl_union_map_align_params(proximity, dim);
2261
2262         if (!domain)
2263                 goto error;
2264
2265         graph.n = isl_union_set_n_set(domain);
2266         if (graph.n == 0)
2267                 goto empty;
2268         if (graph_alloc(ctx, &graph, graph.n,
2269             isl_union_map_n_map(validity) + isl_union_map_n_map(proximity)) < 0)
2270                 goto error;
2271         graph.root = 1;
2272         graph.n = 0;
2273         if (isl_union_set_foreach_set(domain, &extract_node, &graph) < 0)
2274                 goto error;
2275         if (graph_init_table(ctx, &graph) < 0)
2276                 goto error;
2277         graph.n_edge = 0;
2278         if (isl_union_map_foreach_map(validity, &extract_edge, &graph) < 0)
2279                 goto error;
2280         if (graph_init_edge_table(ctx, &graph) < 0)
2281                 goto error;
2282         if (isl_union_map_foreach_map(proximity, &extract_edge, &graph) < 0)
2283                 goto error;
2284
2285         if (compute_schedule(ctx, &graph) < 0)
2286                 goto error;
2287
2288 empty:
2289         sched = extract_schedule(&graph, isl_union_set_get_dim(domain));
2290
2291         graph_free(ctx, &graph);
2292         isl_union_set_free(domain);
2293         isl_union_map_free(validity);
2294         isl_union_map_free(proximity);
2295
2296         return sched;
2297 error:
2298         graph_free(ctx, &graph);
2299         isl_union_set_free(domain);
2300         isl_union_map_free(validity);
2301         isl_union_map_free(proximity);
2302         return NULL;
2303 }
2304
2305 void *isl_schedule_free(__isl_take isl_schedule *sched)
2306 {
2307         int i;
2308         if (!sched)
2309                 return NULL;
2310         for (i = 0; i < sched->n; ++i) {
2311                 isl_map_free(sched->node[i].sched);
2312                 free(sched->node[i].band_end);
2313         }
2314         isl_dim_free(sched->dim);
2315         free(sched);
2316         return NULL;
2317 }
2318
2319 __isl_give isl_union_map *isl_schedule_get_map(__isl_keep isl_schedule *sched)
2320 {
2321         int i;
2322         isl_union_map *umap;
2323
2324         if (!sched)
2325                 return NULL;
2326
2327         umap = isl_union_map_empty(isl_dim_copy(sched->dim));
2328         for (i = 0; i < sched->n; ++i)
2329                 umap = isl_union_map_add_map(umap,
2330                                             isl_map_copy(sched->node[i].sched));
2331
2332         return umap;
2333 }
2334
2335 int isl_schedule_n_band(__isl_keep isl_schedule *sched)
2336 {
2337         return sched ? sched->n_band : 0;
2338 }
2339
2340 /* Construct a mapping that maps each domain to the band in its schedule
2341  * with the specified band index.  Note that bands with the same index
2342  * but for different domains do not need to be related.
2343  */
2344 __isl_give isl_union_map *isl_schedule_get_band(__isl_keep isl_schedule *sched,
2345         unsigned band)
2346 {
2347         int i;
2348         isl_union_map *umap;
2349
2350         if (!sched)
2351                 return NULL;
2352
2353         umap = isl_union_map_empty(isl_dim_copy(sched->dim));
2354         for (i = 0; i < sched->n; ++i) {
2355                 int start, end;
2356                 isl_map *map;
2357
2358                 if (band >= sched->node[i].n_band)
2359                         continue;
2360
2361                 start = band > 0 ? sched->node[i].band_end[band - 1] : 0;
2362                 end = sched->node[i].band_end[band];
2363
2364                 map = isl_map_copy(sched->node[i].sched);
2365
2366                 map = isl_map_project_out(map, isl_dim_out, end,
2367                                           sched->n_total_row - end);
2368                 map = isl_map_project_out(map, isl_dim_out, 0, start);
2369
2370                 umap = isl_union_map_add_map(umap, map);
2371         }
2372
2373         return umap;
2374 }