PaStiX Handbook  6.3.2
order_grids.c
Go to the documentation of this file.
1 /**
2  *
3  * @file order_grids.c
4  *
5  * PaStiX order grids routines
6  *
7  * @copyright 2004-2023 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8  * Univ. Bordeaux. All rights reserved.
9  *
10  * @version 6.3.2
11  * @author Gregoire Pichon
12  * @author Mathieu Faverge
13  * @date 2023-07-21
14  *
15  **/
16 #include "common.h"
17 #include <spm.h>
18 #include "graph/graph.h"
19 #include "order/order_internal.h"
20
21 /**
22  *******************************************************************************
23  *
24  * @ingroup pastix_ordering
25  *
26  * order_grid2D_classic - Orders a separator without any property
27  *
28  *******************************************************************************
29  *
30  * @param[out] peritab
31  * The inverse permutation array
32  *
33  * @param[in] x0
34  * The index of the first vertex in the first direction
35  *
36  * @param[in] xn
37  * The index of the second vertex in the first direction
38  *
39  * @param[in] y0
40  * The index of the first vertex in the second direction
41  *
42  * @param[in] yn
43  * The index of the second vertex in the second direction
44  *
45  * @param[inout] max_number
46  * The larger number to be attributed
47  *
48  * @param[in] ldax
49  * The leading dimension in the first direction
50  *
51  * @param[in] lday
52  * The leading dimension in the second direction
53  *
54  *******************************************************************************/
55 void
57  pastix_int_t x0,
58  pastix_int_t xn,
59  pastix_int_t y0,
60  pastix_int_t yn,
61  pastix_int_t *max_number,
62  pastix_int_t ldax,
63  pastix_int_t lday )
64 {
65  pastix_int_t i, j;
66  pastix_int_t nx = xn-x0;
67  pastix_int_t ny = yn-y0;
68
69  for (i=0; i<nx; i++){
70  for (j=0; j<ny; j++){
71  pastix_int_t index = (x0 + i) * ldax + (y0 + j) * lday;
72  peritab[index] = max_number[0]--;
73  }
74  }
75 }
76
77 /**
78  *******************************************************************************
79  *
80  * @ingroup pastix_ordering
81  *
82  * order_grid3D_classic - Order a 3D Laplacian with quasi optimal ordering:
83  * the current separator is selected as a plan in the smaller direction.
84  *
85  *******************************************************************************
86  *
87  * @param[out] rangtab
88  * The rangtab array that describes the supernodes in the graph.
89  *
90  * @param[out] peritab
91  * The inverse permutation array
92  *
93  * @param[out] cblknbr
94  * The number of supernodes. The internal rangtab and treetab arrays
95  * are of size cblknbr+1
96  *
97  * @param[in] x0
98  * The index of the first vertex in the first direction
99  *
100  * @param[in] xn
101  * The index of the second vertex in the first direction
102  *
103  * @param[in] y0
104  * The index of the first vertex in the second direction
105  *
106  * @param[in] yn
107  * The index of the second vertex in the second direction
108  *
109  * @param[in] z0
110  * The index of the first vertex in the third direction
111  *
112  * @param[in] zn
113  * The index of the second vertex in the third direction
114  *
115  * @param[inout] max_number
116  * The larger number to be attributed
117  *
118  * @param[out] current_rangtab
119  * The index of the current supernode in rangtab
120  *
121  * @param[out] treetab
122  * The treetab array that describes the elimination tree
123  *
124  * @param[in] current_treetab
125  * The index of the current supernode in treetab
126  *
127  * @param[in] ldax
128  * The leading dimension in the first direction
129  *
130  * @param[in] lday
131  * The leading dimension in the second direction
132  *
133  * @param[in] ldaz
134  * The leading dimension in the third direction
135  *
136  *******************************************************************************/
137 void
139  pastix_int_t *peritab,
140  pastix_int_t *cblknbr,
141  pastix_int_t x0,
142  pastix_int_t xn,
143  pastix_int_t y0,
144  pastix_int_t yn,
145  pastix_int_t z0,
146  pastix_int_t zn,
147  pastix_int_t *max_number,
148  pastix_int_t *current_rangtab,
149  pastix_int_t *treetab,
150  pastix_int_t current_treetab,
151  pastix_int_t ldax,
152  pastix_int_t lday,
153  pastix_int_t ldaz )
154 {
155  pastix_int_t dir, i, j;
156  pastix_int_t nx = xn-x0;
157  pastix_int_t ny = yn-y0;
158  pastix_int_t nz = zn-z0;
159  pastix_int_t *peritab_separator;
160
161  /* The subgraph is small enough */
162  if (nx*ny*nz < 15){
163  pastix_int_t k;
164  pastix_int_t current = 0;
165  cblknbr[0] ++;
166  for (i=x0; i<xn; i++){
167  for (j=y0; j<yn; j++){
168  for (k=z0; k<zn; k++){
169  pastix_int_t index = i + ldax * j + ldax*lday * k;
170  peritab[index] = max_number[0] - current;
171  current++;
172  }
173  }
174  }
175
176  treetab[current_rangtab[0]] = current_treetab;
177  rangtab[current_rangtab[0]] = max_number[0];
178  max_number[0] -= current;
179  current_rangtab[0]++;
180  return;
181  }
182
183  cblknbr[0] ++;
184
185  /* In which direction do we cut? 0 for x, 1 for y */
186  dir = 0;
187  if (ny > nx) {
188  dir = 1;
189  }
190  if ((nz > nx) && (nz > ny)) {
191  dir = 2;
192  }
193
194  /* If we cut in direction x */
195  if (dir == 0){
196  treetab[current_rangtab[0]] = current_treetab;
197  rangtab[current_rangtab[0]] = max_number[0];
198  current_rangtab[0]++;
199
200  /* Order separator */
201  peritab_separator = peritab + x0 + nx/2;
202  order_grid2D_classic(peritab_separator,
203  y0, yn, z0, zn,
204  max_number,
205  ldax, ldax * lday);
206
207  /* Order subparts with nested dissection */
208  order_grid3D_classic(rangtab, peritab, cblknbr,
209  x0, x0 + nx/2, y0, yn, z0, zn, max_number,
210  current_rangtab,
211  treetab, current_treetab+1,
212  ldax, lday, ldaz);
213
214  order_grid3D_classic(rangtab, peritab, cblknbr,
215  x0+nx/2+1, xn, y0, yn, z0, zn, max_number,
216  current_rangtab,
217  treetab, current_treetab+1,
218  ldax, lday, ldaz);
219
220  }
221
222  /* If we cut in direction y */
223  else if (dir == 1){
224  treetab[current_rangtab[0]] = current_treetab;
225  rangtab[current_rangtab[0]] = max_number[0];
226  current_rangtab[0]++;
227
228  /* Order separator */
229  peritab_separator = peritab + ldax * (y0 + ny / 2);
230  order_grid2D_classic(peritab_separator,
231  x0, xn, z0, zn,
232  max_number,
233  1, ldax * lday);
234
235  /* Order subparts with nested dissection */
236  order_grid3D_classic(rangtab, peritab, cblknbr,
237  x0, xn, y0, y0+ny/2, z0, zn, max_number,
238  current_rangtab,
239  treetab, current_treetab+1,
240  ldax, lday, ldaz);
241
242  order_grid3D_classic(rangtab, peritab, cblknbr,
243  x0, xn, y0+ny/2+1, yn, z0, zn, max_number,
244  current_rangtab,
245  treetab, current_treetab+1,
246  ldax, lday, ldaz);
247  }
248
249  /* If we cut in direction z */
250  else{
251
252  treetab[current_rangtab[0]] = current_treetab;
253  rangtab[current_rangtab[0]] = max_number[0];
254  current_rangtab[0]++;
255
256  /* Order separator */
257  peritab_separator = peritab + ldax * lday * (z0 + nz/2);
258  order_grid2D_classic(peritab_separator,
259  x0, xn, y0, yn,
260  max_number,
261  1, ldax);
262
263  /* Order subparts with nested dissection */
264  order_grid3D_classic(rangtab, peritab, cblknbr,
265  x0, xn, y0, yn, z0, z0+nz/2, max_number,
266  current_rangtab,
267  treetab, current_treetab+1,
268  ldax, lday, ldaz);
269
270  order_grid3D_classic(rangtab, peritab, cblknbr,
271  x0, xn, y0, yn, z0+nz/2+1, zn, max_number,
272  current_rangtab,
273  treetab, current_treetab+1,
274  ldax, lday, ldaz);
275  }
276 }
277
278 /**
279  *******************************************************************************
280  *
281  * @ingroup pastix_ordering
282  *
283  * orderComputeOptimal - Compute the ordering of a regular 2D or 3D laplacian,
284  * with an optimal strategy.
285  *
286  *******************************************************************************
287  *
288  * @param[out] myorder
289  * Personal ordering for regular laplacians
290  *
291  * @param[in] nx
292  * The number of vertices for the first dimension of the graph
293  *
294  * @param[in] ny
295  * The number of vertices for the second dimension of the graph
296  *
297  * @param[in] nz
298  * The number of vertices for the third dimension of the graph
299  *
300  *******************************************************************************
301  *
302  * @return
303  * \retval PASTIX_SUCCESS on successful exit
304  *
305  *******************************************************************************/
306 int
308  pastix_int_t nx,
309  pastix_int_t ny,
310  pastix_int_t nz )
311 {
312  pastix_order_t *ordemesh = *myorder;
313  pastix_int_t n = nx * ny * nz;
314
315  pastixOrderAlloc(ordemesh, n, n);
316
317  pastix_int_t *rangtab = ordemesh->rangtab;
318  pastix_int_t *permtab = ordemesh->permtab;
319  pastix_int_t *peritab = ordemesh->peritab;
320  pastix_int_t *treetab = ordemesh->treetab;
321
322  pastix_int_t *saved_rangtab, *saved_treetab;
323  pastix_int_t i;
324
325  pastix_int_t current_rangtab = 0;
326
327  /* Graphs for using classical separators */
328  if (nx == ny && ny == nz){
329  pastix_int_t i = 2;
330  while (i != nx && i < nx+1){
331  i = 2*i+1;
332  }
333  if (i != nx){
334  pastix_print_warning( "The given graph size is not correct for optimal manual ordering on 2D regular grid or 3D regular cube. Closer valid sizes are %ld %ld\n",
335  (long)i, (long)(2*i+1));
336  }
337  }
338
339  ordemesh->cblknbr = 0;
340
341  pastix_int_t current_number = n-1;
342  order_grid3D_classic(rangtab, permtab, &ordemesh->cblknbr,
343  0, nx, 0, ny, 0, nz, &current_number, &current_rangtab,
344  treetab, 1,
345  nx, ny, nz);
346
347  for (i=0; i<n; i++){
348  peritab[permtab[i]] = i;
349  }
350
351  saved_rangtab = malloc(n*sizeof(pastix_int_t));
352  memcpy(saved_rangtab, rangtab, n*sizeof(pastix_int_t));
353  saved_treetab = malloc(n*sizeof(pastix_int_t));
354  memcpy(saved_treetab, treetab, n*sizeof(pastix_int_t));
355
356  rangtab[0] = 0;
357  for (i=0; i<ordemesh->cblknbr; i++){
358  rangtab[i+1] = saved_rangtab[ordemesh->cblknbr - i - 1]+1;
359  treetab[i] = saved_treetab[ordemesh->cblknbr - i - 1];
360  }
361  free(saved_rangtab);
362  free(saved_treetab);
363
364  for (i=0; i<ordemesh->cblknbr-1; i++){
365  pastix_int_t j;
366  for (j=i+1; j<ordemesh->cblknbr; j++){
367  if (treetab[j] < treetab[i]){
368  treetab[i] = j;
369  break;
370  }
371  }
372  }
373  treetab[ordemesh->cblknbr-1] = -1;
374
375  ordemesh->rangtab =
376  (pastix_int_t *) memRealloc (rangtab,
377  (ordemesh->cblknbr + 1)*sizeof (pastix_int_t));
378  ordemesh->treetab =
379  (pastix_int_t *) memRealloc (treetab,
380  (ordemesh->cblknbr)*sizeof (pastix_int_t));
381
382  return PASTIX_SUCCESS;
383 }
384
BEGIN_C_DECLS typedef int pastix_int_t
Definition: datatypes.h:51
@ PASTIX_SUCCESS
Definition: api.h:367
pastix_int_t * treetab
Definition: order.h:54
pastix_int_t * permtab
Definition: order.h:51
pastix_int_t * peritab
Definition: order.h:52
pastix_int_t cblknbr
Definition: order.h:50
pastix_int_t * rangtab
Definition: order.h:53
int pastixOrderAlloc(pastix_order_t *ordeptr, pastix_int_t vertnbr, pastix_int_t cblknbr)
Allocate the order structure.
Definition: order.c:55
int pastixOrderGrid(pastix_order_t **myorder, pastix_int_t nx, pastix_int_t ny, pastix_int_t nz)
Definition: order_grids.c:307
Order structure.
Definition: order.h:47
void order_grid2D_classic(pastix_int_t *peritab, pastix_int_t x0, pastix_int_t xn, pastix_int_t y0, pastix_int_t yn, pastix_int_t *max_number, pastix_int_t ldax, pastix_int_t lday)
Definition: order_grids.c:56
void order_grid3D_classic(pastix_int_t *rangtab, pastix_int_t *peritab, pastix_int_t *cblknbr, pastix_int_t x0, pastix_int_t xn, pastix_int_t y0, pastix_int_t yn, pastix_int_t z0, pastix_int_t zn, pastix_int_t *max_number, pastix_int_t *current_rangtab, pastix_int_t *treetab, pastix_int_t current_treetab, pastix_int_t ldax, pastix_int_t lday, pastix_int_t ldaz)
Definition: order_grids.c:138