File: mmetis.c

package info (click to toggle)
parmetis 3.1.1-4
  • links: PTS, VCS
  • area: non-free
  • in suites: jessie, jessie-kfreebsd, wheezy
  • size: 25,620 kB
  • ctags: 2,290
  • sloc: ansic: 27,908; makefile: 220
file content (95 lines) | stat: -rw-r--r-- 3,040 bytes parent folder | download | duplicates (5)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
/*
 * Copyright 1997, Regents of the University of Minnesota
 *
 * mmetis.c
 *
 * This is the entry point of ParMETIS_V3_PartMeshKway
 *
 * Started 10/19/96
 * George
 *
 * $Id: mmetis.c,v 1.8 2003/07/25 04:01:04 karypis Exp $
 *
 */

#include <parmetislib.h>


/***********************************************************************************
* This function is the entry point of the parallel k-way multilevel mesh partitionioner. 
* This function assumes nothing about the mesh distribution.
* It is the general case.
************************************************************************************/
void ParMETIS_V3_PartMeshKway(idxtype *elmdist, idxtype *eptr, idxtype *eind, idxtype *elmwgt, 
                 int *wgtflag, int *numflag, int *ncon, int *ncommonnodes, int *nparts, 
		 float *tpwgts, float *ubvec, int *options, int *edgecut, idxtype *part, 
		 MPI_Comm *comm)
{
  int i, nvtxs, nedges, gnedges, npes, mype;
  idxtype *xadj, *adjncy;
  timer TotalTmr, Mesh2DualTmr, ParMETISTmr;
  CtrlType ctrl;

  /********************************/
  /* Try and take care bad inputs */
  /********************************/
  if (elmdist == NULL || eptr == NULL || eind == NULL || wgtflag == NULL || 
      numflag == NULL || ncon == NULL || ncommonnodes == NULL || nparts == NULL ||
      tpwgts == NULL || ubvec == NULL || options == NULL || edgecut == NULL || 
      part == NULL || comm == NULL) {
    printf("ERROR: One or more required parameters is NULL. Aborting.\n");
    abort();
  }
  if (((*wgtflag)&2) && elmwgt == NULL) {
    printf("ERROR: elmwgt == NULL when vertex weights were specified. Aborting.\n");
    abort();
  }

  
  SetUpCtrl(&ctrl, *nparts, (options[0] == 1 ? options[PMV3_OPTION_DBGLVL] : 0), *comm);
  npes = ctrl.npes;
  mype = ctrl.mype;

  cleartimer(TotalTmr);
  cleartimer(Mesh2DualTmr);
  cleartimer(ParMETISTmr);

  MPI_Barrier(ctrl.comm);
  starttimer(TotalTmr);
  starttimer(Mesh2DualTmr);

  ParMETIS_V3_Mesh2Dual(elmdist, eptr, eind, numflag, ncommonnodes, &xadj, &adjncy, &(ctrl.comm));

  if (ctrl.dbglvl&DBG_INFO) {
    nvtxs = elmdist[mype+1]-elmdist[mype];
    nedges = xadj[nvtxs] + (*numflag == 0 ? 0 : -1);
    rprintf(&ctrl, "Completed Dual Graph -- Nvtxs: %d, Nedges: %d \n", 
            elmdist[npes], GlobalSESum(&ctrl, nedges));
  }

  MPI_Barrier(ctrl.comm);
  stoptimer(Mesh2DualTmr);


  /***********************/
  /* Partition the graph */
  /***********************/
  starttimer(ParMETISTmr);

  ParMETIS_V3_PartKway(elmdist, xadj, adjncy, elmwgt, NULL, wgtflag, numflag, ncon, 
                       nparts, tpwgts, ubvec, options, edgecut, part, &(ctrl.comm));

  MPI_Barrier(ctrl.comm);
  stoptimer(ParMETISTmr);
  stoptimer(TotalTmr);

  IFSET(ctrl.dbglvl, DBG_TIME, PrintTimer(&ctrl, Mesh2DualTmr,	"   Mesh2Dual"));
  IFSET(ctrl.dbglvl, DBG_TIME, PrintTimer(&ctrl, ParMETISTmr,	"    ParMETIS"));
  IFSET(ctrl.dbglvl, DBG_TIME, PrintTimer(&ctrl, TotalTmr,	"       Total"));

  GKfree((void **)&xadj, (void **)&adjncy, LTERM);

  FreeCtrl(&ctrl);

  return;
}