1 /*
2 *
3 * This file is part of MUMPS 4.10.0, built on Tue May 10 12:56:32 UTC 2011
4 *
5 *
6 * This version of MUMPS is provided to you free of charge. It is public
7 * domain, based on public domain software developed during the Esprit IV
8 * European project PARASOL (1996-1999). Since this first public domain
9 * version in 1999, research and developments have been supported by the
10 * following institutions: CERFACS, CNRS, ENS Lyon, INPT(ENSEEIHT)-IRIT,
11 * INRIA, and University of Bordeaux.
12 *
13 * The MUMPS team at the moment of releasing this version includes
14 * Patrick Amestoy, Maurice Bremond, Alfredo Buttari, Abdou Guermouche,
15 * Guillaume Joslin, Jean-Yves L'Excellent, Francois-Henry Rouet, Bora
16 * Ucar and Clement Weisbecker.
17 *
18 * We are also grateful to Emmanuel Agullo, Caroline Bousquet, Indranil
19 * Chowdhury, Philippe Combes, Christophe Daniel, Iain Duff, Vincent Espirat,
20 * Aurelia Fevre, Jacko Koster, Stephane Pralet, Chiara Puglisi, Gregoire
21 * Richard, Tzvetomila Slavova, Miroslav Tuma and Christophe Voemel who
22 * have been contributing to this project.
23 *
24 * Up-to-date copies of the MUMPS package can be obtained
25 * from the Web pages:
26 * http://mumps.enseeiht.fr/ or http://graal.ens-lyon.fr/MUMPS
27 *
28 *
29 * THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY
30 * EXPRESSED OR IMPLIED. ANY USE IS AT YOUR OWN RISK.
31 *
32 *
33 * User documentation of any code that uses this software can
34 * include this complete notice. You can acknowledge (using
35 * references [1] and [2]) the contribution of this package
36 * in any scientific publication dependent upon the use of the
37 * package. You shall use reasonable endeavours to notify
38 * the authors of the package of this publication.
39 *
40 * [1] P. R. Amestoy, I. S. Duff, J. Koster and J.-Y. L'Excellent,
41 * A fully asynchronous multifrontal solver using distributed dynamic
42 * scheduling, SIAM Journal of Matrix Analysis and Applications,
43 * Vol 23, No 1, pp 15-41 (2001).
44 *
45 * [2] P. R. Amestoy and A. Guermouche and J.-Y. L'Excellent and
46 * S. Pralet, Hybrid scheduling for the parallel solution of linear
47 * systems. Parallel Computing Vol 32 (2), pp 136-156 (2006).
48 *
49 */
50 /*
51 * This file contains interfaces to external ordering packages.
52 * At the moment, PORD (J. Schulze) and SCOTCH are interfaced.
53 */
54 #include "mumps_orderings.h"
55 #if defined(pord)
56 /* Interface to PORD */
57 /*int mumps_pord( int, int, int *, int *, int * );
58 #define MUMPS_PORDF \
59 F_SYMBOL(pordf,PORDF)*/
60 void MUMPS_CALL
MUMPS_PORDF(int * nvtx,int * nedges,int * xadj,int * adjncy,int * nv,int * ncmpa)61 MUMPS_PORDF( int *nvtx, int *nedges,
62 int *xadj, int *adjncy,
63 int *nv, int *ncmpa )
64 {
65 *ncmpa = mumps_pord( *nvtx, *nedges, xadj, adjncy, nv );
66 }
67 /* Interface to PORD with weighted graph*/
68 /*int mumps_pord_wnd( int, int, int *, int *, int *, int * );
69 #define MUMPS_PORDF_WND \
70 F_SYMBOL(pordf_wnd,PORDF_WND)*/
71 void MUMPS_CALL
MUMPS_PORDF_WND(int * nvtx,int * nedges,int * xadj,int * adjncy,int * nv,int * ncmpa,int * totw)72 MUMPS_PORDF_WND( int *nvtx, int *nedges,
73 int *xadj, int *adjncy,
74 int *nv, int *ncmpa, int *totw )
75 {
76 *ncmpa = mumps_pord_wnd( *nvtx, *nedges, xadj, adjncy, nv, totw );
77 }
78 /************************************************************
79 mumps_pord is used in ana_aux.F
80 permutation and inverse permutation not set in output,
81 but may be printed in default file: "perm_pord" and "iperm_pord",
82 if associated part uncommneted.
83 But, if uncommetnted a bug occurs in psl_ma41_analysi.F
84 ******************************************************************/
85 /*********************************************************/
mumps_pord(int nvtx,int nedges,int * xadj_pe,int * adjncy,int * nv)86 int mumps_pord
87 (
88 int nvtx,
89 int nedges,
90 int *xadj_pe,
91 int *adjncy,
92 int *nv
93 )
94 {
95 /**********************************
96 Argument Comments:
97 input:
98 -----
99 - nvtx : dimension of the Problem (N)
100 - nedges : number of entries (NZ)
101 - adjncy : non-zeros entries (IW input)
102 input/output:
103 -------------
104 - xadj_pe : pointer through beginning of column non-zeros entries (PTRAR)
105 - on exit, "father array" (PE)
106 ouput:
107 ------
108 - nv : "nfront array" (NV)
109 *************************************/
110 graph_t *G;
111 elimtree_t *T;
112 timings_t cpus[12];
113 options_t options[] = { SPACE_ORDTYPE, SPACE_NODE_SELECTION1,
114 SPACE_NODE_SELECTION2, SPACE_NODE_SELECTION3,
115 SPACE_DOMAIN_SIZE, 0 };
116 int *ncolfactor, *ncolupdate, *parent, *vtx2front;
117 int *first, *link, nfronts, J, K, u, vertex, vertex_root, count;
118 /**************************************************
119 declaration to uncomment if printing ordering
120 ***************************************************
121 FILE *fp1, *fp2;
122 int *perm, *iperm;
123 */
124 /*** decalage des indices couteux dans un premier temps:
125 **** A modifier dans une version ulterieure de MA41GD */
126 for (u = nvtx; u >= 0; u--)
127 {
128 xadj_pe[u] = xadj_pe[u] - 1;
129 }
130 for (K = nedges-1; K >= 0; K--)
131 {
132 adjncy[K] = adjncy[K] - 1;
133 }
134 /* initialization of the graph */
135 mymalloc(G, 1, graph_t);
136 G->xadj = xadj_pe;
137 G->adjncy = adjncy;
138 mymalloc(G->vwght, nvtx, int);
139 G->nvtx = nvtx;
140 G->nedges = nedges;
141 G->type = UNWEIGHTED;
142 G->totvwght = nvtx;
143 for (u = 0; u < nvtx; u++)
144 G->vwght[u] = 1;
145 /* main function of the Ordering */
146 T = SPACE_ordering(G, options, cpus);
147 nfronts = T->nfronts;
148 ncolfactor = T->ncolfactor;
149 ncolupdate = T->ncolupdate;
150 parent = T->parent;
151 /* firstchild = T->firstchild; */
152 vtx2front = T->vtx2front;
153 /* -----------------------------------------------------------
154 store the vertices/columns of a front in a bucket structure
155 ----------------------------------------------------------- */
156 mymalloc(first, nfronts, int);
157 mymalloc(link, nvtx, int);
158 for (J = 0; J < nfronts; J++)
159 first[J] = -1;
160 for (u = nvtx-1; u >= 0; u--)
161 {
162 J = vtx2front[u];
163 link[u] = first[J];
164 first[J] = u;
165 }
166 /* -----------------------------------------------------------
167 fill the two arrays corresponding to the MUMPS tree structure
168 ----------------------------------------------------------- */
169 count = 0;
170 for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
171 {
172 vertex_root = first[K];
173 if (vertex_root == -1)
174 {
175 /* JY: I think this cannot happen */
176 printf(" Internal error in mumps_pord (cf JY), %d\n",K);
177 exit(-1);
178 }
179 /* for the principal column of the supervariable */
180 if (parent[K] == -1)
181 xadj_pe[vertex_root] = 0; /* root of the tree */
182 else
183 xadj_pe[vertex_root] = - (first[parent[K]]+1);
184 nv[vertex_root] = ncolfactor[K] + ncolupdate[K];
185 count++;
186 for (vertex = link[vertex_root]; vertex != -1; vertex = link[vertex])
187 /* for the secondary columns of the supervariable */
188 {
189 xadj_pe[vertex] = - (vertex_root+1);
190 nv[vertex] = 0;
191 count++;
192 }
193 }
194 /* ----------------------
195 free memory and return
196 ---------------------- */
197 free(first); free(link);
198 free(G->vwght);
199 free(G);
200 freeElimTree(T);
201 return (0);
202 }
203 /*********************************************************/
mumps_pord_wnd(int nvtx,int nedges,int * xadj_pe,int * adjncy,int * nv,int * totw)204 int mumps_pord_wnd
205 (
206 int nvtx,
207 int nedges,
208 int *xadj_pe,
209 int *adjncy,
210 int *nv,
211 int *totw
212 )
213 {
214 /**********************************
215 Argument Comments:
216 input:
217 -----
218 - nvtx : dimension of the Problem (N)
219 - nedges : number of entries (NZ)
220 - adjncy : non-zeros entries (IW input)
221 - totw : sum of the weigth of the vertices
222 input/output:
223 -------------
224 - xadj_pe : pointer through beginning of column non-zeros entries (PTRAR)
225 - on exit, "father array" (PE)
226 ouput:
227 ------
228 - nv : weight of the vertices
229 - on exit "nfront array" (NV)
230 *************************************/
231 graph_t *G;
232 elimtree_t *T;
233 timings_t cpus[12];
234 options_t options[] = { SPACE_ORDTYPE, SPACE_NODE_SELECTION1,
235 SPACE_NODE_SELECTION2, SPACE_NODE_SELECTION3,
236 SPACE_DOMAIN_SIZE, 0 };
237 int *ncolfactor, *ncolupdate, *parent, *vtx2front;
238 int *first, *link, nfronts, J, K, u, vertex, vertex_root, count;
239 /**************************************************
240 declaration to uncomment if printing ordering
241 ***************************************************
242 FILE *fp1, *fp2;
243 int *perm, *iperm;
244 */
245 /*** decalage des indices couteux dans un premier temps:
246 **** A modifier dans une version ulterieure de MA41GD */
247 for (u = nvtx; u >= 0; u--)
248 {
249 xadj_pe[u] = xadj_pe[u] - 1;
250 }
251 for (K = nedges-1; K >= 0; K--)
252 {
253 adjncy[K] = adjncy[K] - 1;
254 }
255 /* initialization of the graph */
256 mymalloc(G, 1, graph_t);
257 G->xadj = xadj_pe;
258 G->adjncy= adjncy;
259 mymalloc(G->vwght, nvtx, int);
260 G->nvtx = nvtx;
261 G->nedges = nedges;
262 G->type = WEIGHTED;
263 G->totvwght = (*totw);
264 for (u = 0; u < nvtx; u++)
265 G->vwght[u] = nv[u];
266 /* main function of the Ordering */
267 T = SPACE_ordering(G, options, cpus);
268 nfronts = T->nfronts;
269 ncolfactor = T->ncolfactor;
270 ncolupdate = T->ncolupdate;
271 parent = T->parent;
272 /* firstchild = T->firstchild; */
273 vtx2front = T->vtx2front;
274 /* -----------------------------------------------------------
275 store the vertices/columns of a front in a bucket structure
276 ----------------------------------------------------------- */
277 mymalloc(first, nfronts, int);
278 mymalloc(link, nvtx, int);
279 for (J = 0; J < nfronts; J++)
280 first[J] = -1;
281 for (u = nvtx-1; u >= 0; u--)
282 {
283 J = vtx2front[u];
284 link[u] = first[J];
285 first[J] = u;
286 }
287 /* -----------------------------------------------------------
288 fill the two arrays corresponding to the MUMPS tree structure
289 ----------------------------------------------------------- */
290 count = 0;
291 for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
292 {
293 vertex_root = first[K];
294 if (vertex_root == -1)
295 {
296 /* JY: I think this cannot happen */
297 printf(" Internal error in mumps_pord (cf JY), %d\n",K);
298 exit(-1);
299 }
300 /* for the principal column of the supervariable */
301 if (parent[K] == -1)
302 xadj_pe[vertex_root] = 0; /* root of the tree */
303 else
304 xadj_pe[vertex_root] = - (first[parent[K]]+1);
305 nv[vertex_root] = ncolfactor[K] + ncolupdate[K];
306 count++;
307 for (vertex = link[vertex_root]; vertex != -1; vertex = link[vertex])
308 /* for the secondary columns of the supervariable */
309 {
310 xadj_pe[vertex] = - (vertex_root+1);
311 nv[vertex] = 0;
312 count++;
313 }
314 }
315 /* ----------------------
316 free memory and return
317 ---------------------- */
318 free(first); free(link);
319 free(G->vwght);
320 free(G);
321 freeElimTree(T);
322 return (0);
323 }
324 #endif /* pord */
325 /************************************************************/
326 #if defined(scotch) || defined(ptscotch)
327 /*int esmumps( const int n, const int iwlen, int * const pe, const int pfree,
328 int * const len, int * const iw, int * const nv, int * const elen,
329 int * const last);*/
330 /* Fortran interface to SCOTCH */
331 /*#define MUMPS_SCOTCH \
332 F_SYMBOL(scotch,SCOTCH)*/
333 void MUMPS_CALL
MUMPS_SCOTCH(const int * const n,const int * const iwlen,int * const petab,const int * const pfree,int * const lentab,int * const iwtab,int * const nvtab,int * const elentab,int * const lasttab,int * const ncmpa)334 MUMPS_SCOTCH( const int * const n,
335 const int * const iwlen,
336 int * const petab,
337 const int * const pfree,
338 int * const lentab,
339 int * const iwtab,
340 int * const nvtab,
341 int * const elentab,
342 int * const lasttab,
343 int * const ncmpa )
344 {
345 *ncmpa = esmumps( *n, *iwlen, petab, *pfree,
346 lentab, iwtab, nvtab, elentab, lasttab );
347 }
348 #endif /* scotch */
349 #if defined(ptscotch)
350 /*#include "mpi.h"
351 #include <stdio.h>
352 #include "ptscotch.h"
353 int mumps_dgraphinit( SCOTCH_Dgraph *, MPI_Fint *, MPI_Fint *);
354 #define MUMPS_DGRAPHINIT \
355 F_SYMBOL(dgraphinit,DGRAPHINIT)*/
356 void MUMPS_CALL
MUMPS_DGRAPHINIT(SCOTCH_Dgraph * graphptr,MPI_Fint * comm,MPI_Fint * ierr)357 MUMPS_DGRAPHINIT(SCOTCH_Dgraph *graphptr, MPI_Fint *comm, MPI_Fint *ierr)
358 {
359 MPI_Comm int_comm;
360 int_comm = MPI_Comm_f2c(*comm);
361 *ierr = SCOTCH_dgraphInit(graphptr, int_comm);
362 return;
363 }
364 #endif
365 #if defined(parmetis)
366 void MUMPS_CALL
MUMPS_PARMETIS(int * first,int * vertloctab,int * edgeloctab,int * numflag,int * options,int * order,int * sizes,int * comm)367 MUMPS_PARMETIS(int *first, int *vertloctab,
368 int *edgeloctab, int *numflag,
369 int *options, int *order,
370 int *sizes, int *comm)
371 {
372 MPI_Comm int_comm;
373 int_comm = MPI_Comm_f2c(*comm);
374 ParMETIS_V3_NodeND(first, vertloctab, edgeloctab, numflag, options, order, sizes, &int_comm);
375 return;
376 }
377 #endif
378