Charm++: Charm Source Code Documentation

00001 /*
00002  * Copyright 1997, Regents of the University of Minnesota
00003  *
00004  * mesh.c
00005  *
00006  * This file contains routines for converting 3D and 4D finite element
00007  * meshes into dual or nodal graphs
00008  *
00009  * Started 8/18/97
00010  * George
00011  *
00012  * $Id$
00013  *
00014  */
00015 
00016 #include <metis.h>
00017 
00018 /*****************************************************************************
00019 * This function creates a graph corresponding to the dual of a finite element
00020 * mesh. At this point the supported elements are triangles, tetrahedrons, and
00021 * bricks.
00022 ******************************************************************************/
00023 void METIS_MeshToDual(int *ne, int *nn, idxtype *elmnts, int *etype, int *numflag, 
00024                       idxtype *dxadj, idxtype *dadjncy)
00025 {
00026   int esizes[] = {-1, 3, 4, 8, 4};
00027 
00028   if (*numflag == 1)
00029     ChangeMesh2CNumbering((*ne)*esizes[*etype], elmnts);
00030 
00031   GENDUALMETIS(*ne, *nn, *etype, elmnts, dxadj, dadjncy);
00032 
00033   if (*numflag == 1)
00034     ChangeMesh2FNumbering((*ne)*esizes[*etype], elmnts, *ne, dxadj, dadjncy);
00035 }
00036 
00037 
00038 /*****************************************************************************
00039 * This function creates a graph corresponding to the finite element mesh. 
00040 * At this point the supported elements are triangles, tetrahedrons.
00041 ******************************************************************************/
00042 void METIS_MeshToNodal(int *ne, int *nn, idxtype *elmnts, int *etype, int *numflag, 
00043                        idxtype *dxadj, idxtype *dadjncy)
00044 {
00045   int esizes[] = {-1, 3, 4, 8, 4};
00046 
00047   if (*numflag == 1)
00048     ChangeMesh2CNumbering((*ne)*esizes[*etype], elmnts);
00049 
00050   switch (*etype) {
00051     case 1:
00052       TRINODALMETIS(*ne, *nn, elmnts, dxadj, dadjncy);
00053       break;
00054     case 2:
00055       TETNODALMETIS(*ne, *nn, elmnts, dxadj, dadjncy);
00056       break;
00057     case 3:
00058       HEXNODALMETIS(*ne, *nn, elmnts, dxadj, dadjncy);
00059       break;
00060     case 4:
00061       QUADNODALMETIS(*ne, *nn, elmnts, dxadj, dadjncy);
00062       break;
00063   }
00064 
00065   if (*numflag == 1)
00066     ChangeMesh2FNumbering((*ne)*esizes[*etype], elmnts, *nn, dxadj, dadjncy);
00067 }
00068 
00069 
00070 
00071 /*****************************************************************************
00072 * This function creates the dual of a finite element mesh
00073 ******************************************************************************/
00074 void GENDUALMETIS(int nelmnts, int nvtxs, int etype, idxtype *elmnts, idxtype *dxadj, 
00075                   idxtype *dadjncy)
00076 {
00077    int i, j, jj, k, kk, kkk, l, m, n, nedges, mask;
00078    idxtype *nptr, *nind;
00079    idxtype *mark, ind[200], wgt[200];
00080    int esize, esizes[] = {-1, 3, 4, 8, 4},
00081        mgcnum, mgcnums[] = {-1, 2, 3, 4, 2};
00082 
00083    mask = (1<<11)-1;
00084    mark = idxsmalloc(mask+1, -1, "GENDUALMETIS: mark");
00085 
00086    /* Get the element size and magic number for the particular element */
00087    esize = esizes[etype];
00088    mgcnum = mgcnums[etype];
00089 
00090    /* Construct the node-element list first */
00091    nptr = idxsmalloc(nvtxs+1, 0, "GENDUALMETIS: nptr");
00092    for (j=esize*nelmnts, i=0; i<j; i++) 
00093      nptr[elmnts[i]]++;
00094    MAKECSR(i, nvtxs, nptr);
00095 
00096    nind = idxmalloc(nptr[nvtxs], "GENDUALMETIS: nind");
00097    for (k=i=0; i<nelmnts; i++) {
00098      for (j=0; j<esize; j++, k++) 
00099        nind[nptr[elmnts[k]]++] = i;
00100    }
00101    for (i=nvtxs; i>0; i--)
00102      nptr[i] = nptr[i-1];
00103    nptr[0] = 0;
00104 
00105    for (i=0; i<nelmnts; i++) 
00106      dxadj[i] = esize*i;
00107 
00108    for (i=0; i<nelmnts; i++) {
00109      for (m=j=0; j<esize; j++) {
00110        n = elmnts[esize*i+j];
00111        for (k=nptr[n+1]-1; k>=nptr[n]; k--) {
00112          if ((kk = nind[k]) <= i)
00113            break;
00114 
00115          kkk = kk&mask;
00116          if ((l = mark[kkk]) == -1) {
00117            ind[m] = kk;
00118            wgt[m] = 1;
00119            mark[kkk] = m++;
00120          }
00121          else if (ind[l] == kk) {
00122            wgt[l]++;
00123          }
00124          else {
00125            for (jj=0; jj<m; jj++) {
00126              if (ind[jj] == kk) {
00127                wgt[jj]++;
00128                break;
00129              }
00130            }
00131            if (jj == m) {
00132              ind[m] = kk;
00133              wgt[m++] = 1;
00134            }
00135          }
00136        }
00137      }
00138      for (j=0; j<m; j++) {
00139        if (wgt[j] == mgcnum) {
00140          k = ind[j];
00141          dadjncy[dxadj[i]++] = k;
00142          dadjncy[dxadj[k]++] = i;
00143        }
00144        mark[ind[j]&mask] = -1;
00145      }
00146    }
00147 
00148    /* Go and consolidate the dxadj and dadjncy */
00149    for (j=i=0; i<nelmnts; i++) {
00150      for (k=esize*i; k<dxadj[i]; k++, j++)
00151        dadjncy[j] = dadjncy[k];
00152      dxadj[i] = j;
00153    }
00154    for (i=nelmnts; i>0; i--)
00155      dxadj[i] = dxadj[i-1];
00156    dxadj[0] = 0;
00157 
00158    free(mark);
00159    free(nptr);
00160    free(nind);
00161 
00162 }
00163 
00164 
00165 
00166 
00167 /*****************************************************************************
00168 * This function creates the nodal graph of a finite element mesh
00169 ******************************************************************************/
00170 void TRINODALMETIS(int nelmnts, int nvtxs, idxtype *elmnts, idxtype *dxadj, idxtype *dadjncy)
00171 {
00172    int i, j, jj, k, kk, kkk, l, m, n, nedges;
00173    idxtype *nptr, *nind;
00174    idxtype *mark;
00175 
00176    /* Construct the node-element list first */
00177    nptr = idxsmalloc(nvtxs+1, 0, "TRINODALMETIS: nptr");
00178    for (j=3*nelmnts, i=0; i<j; i++) 
00179      nptr[elmnts[i]]++;
00180    MAKECSR(i, nvtxs, nptr);
00181 
00182    nind = idxmalloc(nptr[nvtxs], "TRINODALMETIS: nind");
00183    for (k=i=0; i<nelmnts; i++) {
00184      for (j=0; j<3; j++, k++) 
00185        nind[nptr[elmnts[k]]++] = i;
00186    }
00187    for (i=nvtxs; i>0; i--)
00188      nptr[i] = nptr[i-1];
00189    nptr[0] = 0;
00190 
00191 
00192    mark = idxsmalloc(nvtxs, -1, "TRINODALMETIS: mark");
00193 
00194    nedges = dxadj[0] = 0;
00195    for (i=0; i<nvtxs; i++) {
00196      mark[i] = i;
00197      for (j=nptr[i]; j<nptr[i+1]; j++) {
00198        for (jj=3*nind[j], k=0; k<3; k++, jj++) {
00199          kk = elmnts[jj];
00200          if (mark[kk] != i) {
00201            mark[kk] = i;
00202            dadjncy[nedges++] = kk;
00203          }
00204        }
00205      }
00206      dxadj[i+1] = nedges;
00207    }
00208 
00209    free(mark);
00210    free(nptr);
00211    free(nind);
00212 
00213 }
00214 
00215 
00216 /*****************************************************************************
00217 * This function creates the nodal graph of a finite element mesh
00218 ******************************************************************************/
00219 void TETNODALMETIS(int nelmnts, int nvtxs, idxtype *elmnts, idxtype *dxadj, idxtype *dadjncy)
00220 {
00221    int i, j, jj, k, kk, kkk, l, m, n, nedges;
00222    idxtype *nptr, *nind;
00223    idxtype *mark;
00224 
00225    /* Construct the node-element list first */
00226    nptr = idxsmalloc(nvtxs+1, 0, "TETNODALMETIS: nptr");
00227    for (j=4*nelmnts, i=0; i<j; i++) 
00228      nptr[elmnts[i]]++;
00229    MAKECSR(i, nvtxs, nptr);
00230 
00231    nind = idxmalloc(nptr[nvtxs], "TETNODALMETIS: nind");
00232    for (k=i=0; i<nelmnts; i++) {
00233      for (j=0; j<4; j++, k++) 
00234        nind[nptr[elmnts[k]]++] = i;
00235    }
00236    for (i=nvtxs; i>0; i--)
00237      nptr[i] = nptr[i-1];
00238    nptr[0] = 0;
00239 
00240 
00241    mark = idxsmalloc(nvtxs, -1, "TETNODALMETIS: mark");
00242 
00243    nedges = dxadj[0] = 0;
00244    for (i=0; i<nvtxs; i++) {
00245      mark[i] = i;
00246      for (j=nptr[i]; j<nptr[i+1]; j++) {
00247        for (jj=4*nind[j], k=0; k<4; k++, jj++) {
00248          kk = elmnts[jj];
00249          if (mark[kk] != i) {
00250            mark[kk] = i;
00251            dadjncy[nedges++] = kk;
00252          }
00253        }
00254      }
00255      dxadj[i+1] = nedges;
00256    }
00257 
00258    free(mark);
00259    free(nptr);
00260    free(nind);
00261 
00262 }
00263 
00264 
00265 /*****************************************************************************
00266 * This function creates the nodal graph of a finite element mesh
00267 ******************************************************************************/
00268 void HEXNODALMETIS(int nelmnts, int nvtxs, idxtype *elmnts, idxtype *dxadj, idxtype *dadjncy)
00269 {
00270    int i, j, jj, k, kk, kkk, l, m, n, nedges;
00271    idxtype *nptr, *nind;
00272    idxtype *mark;
00273    int table[8][3] = {1, 3, 4,
00274                       0, 2, 5,
00275                       1, 3, 6,
00276                       0, 2, 7,
00277                       0, 5, 7,
00278                       1, 4, 6,
00279                       2, 5, 7,
00280                       3, 4, 6};
00281 
00282    /* Construct the node-element list first */
00283    nptr = idxsmalloc(nvtxs+1, 0, "HEXNODALMETIS: nptr");
00284    for (j=8*nelmnts, i=0; i<j; i++) 
00285      nptr[elmnts[i]]++;
00286    MAKECSR(i, nvtxs, nptr);
00287 
00288    nind = idxmalloc(nptr[nvtxs], "HEXNODALMETIS: nind");
00289    for (k=i=0; i<nelmnts; i++) {
00290      for (j=0; j<8; j++, k++) 
00291        nind[nptr[elmnts[k]]++] = i;
00292    }
00293    for (i=nvtxs; i>0; i--)
00294      nptr[i] = nptr[i-1];
00295    nptr[0] = 0;
00296 
00297 
00298    mark = idxsmalloc(nvtxs, -1, "HEXNODALMETIS: mark");
00299 
00300    nedges = dxadj[0] = 0;
00301    for (i=0; i<nvtxs; i++) {
00302      mark[i] = i;
00303      for (j=nptr[i]; j<nptr[i+1]; j++) {
00304        jj=8*nind[j];
00305        for (k=0; k<8; k++) {
00306          if (elmnts[jj+k] == i)
00307            break;
00308        }
00309        ASSERT(k != 8);
00310 
00311        /* You found the index, now go and put the 3 neighbors */
00312        kk = elmnts[jj+table[k][0]];
00313        if (mark[kk] != i) {
00314          mark[kk] = i;
00315          dadjncy[nedges++] = kk;
00316        }
00317        kk = elmnts[jj+table[k][1]];
00318        if (mark[kk] != i) {
00319          mark[kk] = i;
00320          dadjncy[nedges++] = kk;
00321        }
00322        kk = elmnts[jj+table[k][2]];
00323        if (mark[kk] != i) {
00324          mark[kk] = i;
00325          dadjncy[nedges++] = kk;
00326        }
00327      }
00328      dxadj[i+1] = nedges;
00329    }
00330 
00331    free(mark);
00332    free(nptr);
00333    free(nind);
00334 
00335 }
00336 
00337 
00338 /*****************************************************************************
00339 * This function creates the nodal graph of a finite element mesh
00340 ******************************************************************************/
00341 void QUADNODALMETIS(int nelmnts, int nvtxs, idxtype *elmnts, idxtype *dxadj, idxtype *dadjncy)
00342 {
00343    int i, j, jj, k, kk, kkk, l, m, n, nedges;
00344    idxtype *nptr, *nind;
00345    idxtype *mark;
00346    int table[4][2] = {1, 3, 
00347                       0, 2,
00348                       1, 3, 
00349                       0, 2}; 
00350 
00351    /* Construct the node-element list first */
00352    nptr = idxsmalloc(nvtxs+1, 0, "QUADNODALMETIS: nptr");
00353    for (j=4*nelmnts, i=0; i<j; i++) 
00354      nptr[elmnts[i]]++;
00355    MAKECSR(i, nvtxs, nptr);
00356 
00357    nind = idxmalloc(nptr[nvtxs], "QUADNODALMETIS: nind");
00358    for (k=i=0; i<nelmnts; i++) {
00359      for (j=0; j<4; j++, k++) 
00360        nind[nptr[elmnts[k]]++] = i;
00361    }
00362    for (i=nvtxs; i>0; i--)
00363      nptr[i] = nptr[i-1];
00364    nptr[0] = 0;
00365 
00366 
00367    mark = idxsmalloc(nvtxs, -1, "QUADNODALMETIS: mark");
00368 
00369    nedges = dxadj[0] = 0;
00370    for (i=0; i<nvtxs; i++) {
00371      mark[i] = i;
00372      for (j=nptr[i]; j<nptr[i+1]; j++) {
00373        jj=4*nind[j];
00374        for (k=0; k<4; k++) {
00375          if (elmnts[jj+k] == i)
00376            break;
00377        }
00378        ASSERT(k != 4);
00379 
00380        /* You found the index, now go and put the 2 neighbors */
00381        kk = elmnts[jj+table[k][0]];
00382        if (mark[kk] != i) {
00383          mark[kk] = i;
00384          dadjncy[nedges++] = kk;
00385        }
00386        kk = elmnts[jj+table[k][1]];
00387        if (mark[kk] != i) {
00388          mark[kk] = i;
00389          dadjncy[nedges++] = kk;
00390        }
00391      }
00392      dxadj[i+1] = nedges;
00393    }
00394 
00395    free(mark);
00396    free(nptr);
00397    free(nind);
00398 
00399 }
libs/ck-libs/parmetis/METISLib/mesh.c
