You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
412 lines
12 KiB
412 lines
12 KiB
/*
|
|
* Copyright 1997, Regents of the University of Minnesota
|
|
*
|
|
* mesh.c
|
|
*
|
|
* This file contains routines for converting 3D and 4D finite element
|
|
* meshes into dual or nodal graphs
|
|
*
|
|
* Started 8/18/97
|
|
* George
|
|
*
|
|
* $Id: mesh.c 13804 2013-03-04 23:49:08Z karypis $
|
|
*
|
|
*/
|
|
|
|
#include "metislib.h"
|
|
|
|
|
|
/*****************************************************************************/
|
|
/*! This function creates a graph corresponding to the dual of a finite element
|
|
mesh.
|
|
|
|
\param ne is the number of elements in the mesh.
|
|
\param nn is the number of nodes in the mesh.
|
|
\param eptr is an array of size ne+1 used to mark the start and end
|
|
locations in the nind array.
|
|
\param eind is an array that stores for each element the set of node IDs
|
|
(indices) that it is made off. The length of this array is equal
|
|
to the total number of nodes over all the mesh elements.
|
|
\param ncommon is the minimum number of nodes that two elements must share
|
|
in order to be connected via an edge in the dual graph.
|
|
\param numflag is either 0 or 1 indicating if the numbering of the nodes
|
|
starts from 0 or 1, respectively. The same numbering is used for the
|
|
returned graph as well.
|
|
\param r_xadj indicates where the adjacency list of each vertex is stored
|
|
in r_adjncy. The memory for this array is allocated by this routine.
|
|
It can be freed by calling METIS_free().
|
|
\param r_adjncy stores the adjacency list of each vertex in the generated
|
|
dual graph. The memory for this array is allocated by this routine.
|
|
It can be freed by calling METIS_free().
|
|
|
|
*/
|
|
/*****************************************************************************/
|
|
int METIS_MeshToDual(idx_t *ne, idx_t *nn, idx_t *eptr, idx_t *eind,
|
|
idx_t *ncommon, idx_t *numflag, idx_t **r_xadj, idx_t **r_adjncy)
|
|
{
|
|
int sigrval=0, renumber=0;
|
|
|
|
/* set up malloc cleaning code and signal catchers */
|
|
if (!gk_malloc_init())
|
|
return METIS_ERROR_MEMORY;
|
|
|
|
gk_sigtrap();
|
|
|
|
if ((sigrval = gk_sigcatch()) != 0)
|
|
goto SIGTHROW;
|
|
|
|
|
|
/* renumber the mesh */
|
|
if (*numflag == 1) {
|
|
ChangeMesh2CNumbering(*ne, eptr, eind);
|
|
renumber = 1;
|
|
}
|
|
|
|
/* create dual graph */
|
|
*r_xadj = *r_adjncy = NULL;
|
|
CreateGraphDual(*ne, *nn, eptr, eind, *ncommon, r_xadj, r_adjncy);
|
|
|
|
|
|
SIGTHROW:
|
|
if (renumber)
|
|
ChangeMesh2FNumbering(*ne, eptr, eind, *ne, *r_xadj, *r_adjncy);
|
|
|
|
gk_siguntrap();
|
|
gk_malloc_cleanup(0);
|
|
|
|
if (sigrval != 0) {
|
|
if (*r_xadj != NULL)
|
|
free(*r_xadj);
|
|
if (*r_adjncy != NULL)
|
|
free(*r_adjncy);
|
|
*r_xadj = *r_adjncy = NULL;
|
|
}
|
|
|
|
return metis_rcode(sigrval);
|
|
}
|
|
|
|
|
|
/*****************************************************************************/
|
|
/*! This function creates a graph corresponding to (almost) the nodal of a
|
|
finite element mesh. In the nodal graph, each node is connected to the
|
|
nodes corresponding to the union of nodes present in all the elements
|
|
in which that node belongs.
|
|
|
|
\param ne is the number of elements in the mesh.
|
|
\param nn is the number of nodes in the mesh.
|
|
\param eptr is an array of size ne+1 used to mark the start and end
|
|
locations in the nind array.
|
|
\param eind is an array that stores for each element the set of node IDs
|
|
(indices) that it is made off. The length of this array is equal
|
|
to the total number of nodes over all the mesh elements.
|
|
\param numflag is either 0 or 1 indicating if the numbering of the nodes
|
|
starts from 0 or 1, respectively. The same numbering is used for the
|
|
returned graph as well.
|
|
\param r_xadj indicates where the adjacency list of each vertex is stored
|
|
in r_adjncy. The memory for this array is allocated by this routine.
|
|
It can be freed by calling METIS_free().
|
|
\param r_adjncy stores the adjacency list of each vertex in the generated
|
|
dual graph. The memory for this array is allocated by this routine.
|
|
It can be freed by calling METIS_free().
|
|
|
|
*/
|
|
/*****************************************************************************/
|
|
int METIS_MeshToNodal(idx_t *ne, idx_t *nn, idx_t *eptr, idx_t *eind,
|
|
idx_t *numflag, idx_t **r_xadj, idx_t **r_adjncy)
|
|
{
|
|
int sigrval=0, renumber=0;
|
|
|
|
/* set up malloc cleaning code and signal catchers */
|
|
if (!gk_malloc_init())
|
|
return METIS_ERROR_MEMORY;
|
|
|
|
gk_sigtrap();
|
|
|
|
if ((sigrval = gk_sigcatch()) != 0)
|
|
goto SIGTHROW;
|
|
|
|
|
|
/* renumber the mesh */
|
|
if (*numflag == 1) {
|
|
ChangeMesh2CNumbering(*ne, eptr, eind);
|
|
renumber = 1;
|
|
}
|
|
|
|
/* create nodal graph */
|
|
*r_xadj = *r_adjncy = NULL;
|
|
CreateGraphNodal(*ne, *nn, eptr, eind, r_xadj, r_adjncy);
|
|
|
|
|
|
SIGTHROW:
|
|
if (renumber)
|
|
ChangeMesh2FNumbering(*ne, eptr, eind, *nn, *r_xadj, *r_adjncy);
|
|
|
|
gk_siguntrap();
|
|
gk_malloc_cleanup(0);
|
|
|
|
if (sigrval != 0) {
|
|
if (*r_xadj != NULL)
|
|
free(*r_xadj);
|
|
if (*r_adjncy != NULL)
|
|
free(*r_adjncy);
|
|
*r_xadj = *r_adjncy = NULL;
|
|
}
|
|
|
|
return metis_rcode(sigrval);
|
|
}
|
|
|
|
|
|
/*****************************************************************************/
|
|
/*! This function creates the dual of a finite element mesh */
|
|
/*****************************************************************************/
|
|
void CreateGraphDual(idx_t ne, idx_t nn, idx_t *eptr, idx_t *eind, idx_t ncommon,
|
|
idx_t **r_xadj, idx_t **r_adjncy)
|
|
{
|
|
idx_t i, j, nnbrs;
|
|
idx_t *nptr, *nind;
|
|
idx_t *xadj, *adjncy;
|
|
idx_t *marker, *nbrs;
|
|
|
|
if (ncommon < 1) {
|
|
printf(" Increased ncommon to 1, as it was initially %"PRIDX"\n", ncommon);
|
|
ncommon = 1;
|
|
}
|
|
|
|
/* construct the node-element list first */
|
|
nptr = ismalloc(nn+1, 0, "CreateGraphDual: nptr");
|
|
nind = imalloc(eptr[ne], "CreateGraphDual: nind");
|
|
|
|
for (i=0; i<ne; i++) {
|
|
for (j=eptr[i]; j<eptr[i+1]; j++)
|
|
nptr[eind[j]]++;
|
|
}
|
|
MAKECSR(i, nn, nptr);
|
|
|
|
for (i=0; i<ne; i++) {
|
|
for (j=eptr[i]; j<eptr[i+1]; j++)
|
|
nind[nptr[eind[j]]++] = i;
|
|
}
|
|
SHIFTCSR(i, nn, nptr);
|
|
|
|
|
|
/* Allocate memory for xadj, since you know its size.
|
|
These are done using standard malloc as they are returned
|
|
to the calling function */
|
|
if ((xadj = (idx_t *)malloc((ne+1)*sizeof(idx_t))) == NULL)
|
|
gk_errexit(SIGMEM, "***Failed to allocate memory for xadj.\n");
|
|
*r_xadj = xadj;
|
|
iset(ne+1, 0, xadj);
|
|
|
|
/* allocate memory for working arrays used by FindCommonElements */
|
|
marker = ismalloc(ne, 0, "CreateGraphDual: marker");
|
|
nbrs = imalloc(ne, "CreateGraphDual: nbrs");
|
|
|
|
for (i=0; i<ne; i++) {
|
|
xadj[i] = FindCommonElements(i, eptr[i+1]-eptr[i], eind+eptr[i], nptr,
|
|
nind, eptr, ncommon, marker, nbrs);
|
|
}
|
|
MAKECSR(i, ne, xadj);
|
|
|
|
/* Allocate memory for adjncy, since you now know its size.
|
|
These are done using standard malloc as they are returned
|
|
to the calling function */
|
|
if ((adjncy = (idx_t *)malloc(xadj[ne]*sizeof(idx_t))) == NULL) {
|
|
free(xadj);
|
|
*r_xadj = NULL;
|
|
gk_errexit(SIGMEM, "***Failed to allocate memory for adjncy.\n");
|
|
}
|
|
*r_adjncy = adjncy;
|
|
|
|
for (i=0; i<ne; i++) {
|
|
nnbrs = FindCommonElements(i, eptr[i+1]-eptr[i], eind+eptr[i], nptr,
|
|
nind, eptr, ncommon, marker, nbrs);
|
|
for (j=0; j<nnbrs; j++)
|
|
adjncy[xadj[i]++] = nbrs[j];
|
|
}
|
|
SHIFTCSR(i, ne, xadj);
|
|
|
|
gk_free((void **)&nptr, &nind, &marker, &nbrs, LTERM);
|
|
}
|
|
|
|
|
|
/*****************************************************************************/
|
|
/*! This function finds all elements that share at least ncommon nodes with
|
|
the ``query'' element.
|
|
*/
|
|
/*****************************************************************************/
|
|
idx_t FindCommonElements(idx_t qid, idx_t elen, idx_t *eind, idx_t *nptr,
|
|
idx_t *nind, idx_t *eptr, idx_t ncommon, idx_t *marker, idx_t *nbrs)
|
|
{
|
|
idx_t i, ii, j, jj, k, l, overlap;
|
|
|
|
/* find all elements that share at least one node with qid */
|
|
for (k=0, i=0; i<elen; i++) {
|
|
j = eind[i];
|
|
for (ii=nptr[j]; ii<nptr[j+1]; ii++) {
|
|
jj = nind[ii];
|
|
|
|
if (marker[jj] == 0)
|
|
nbrs[k++] = jj;
|
|
marker[jj]++;
|
|
}
|
|
}
|
|
|
|
/* put qid into the neighbor list (in case it is not there) so that it
|
|
will be removed in the next step */
|
|
if (marker[qid] == 0)
|
|
nbrs[k++] = qid;
|
|
marker[qid] = 0;
|
|
|
|
/* compact the list to contain only those with at least ncommon nodes */
|
|
for (j=0, i=0; i<k; i++) {
|
|
overlap = marker[l = nbrs[i]];
|
|
if (overlap >= ncommon ||
|
|
overlap >= elen-1 ||
|
|
overlap >= eptr[l+1]-eptr[l]-1)
|
|
nbrs[j++] = l;
|
|
marker[l] = 0;
|
|
}
|
|
|
|
return j;
|
|
}
|
|
|
|
|
|
/*****************************************************************************/
|
|
/*! This function creates the (almost) nodal of a finite element mesh */
|
|
/*****************************************************************************/
|
|
void CreateGraphNodal(idx_t ne, idx_t nn, idx_t *eptr, idx_t *eind,
|
|
idx_t **r_xadj, idx_t **r_adjncy)
|
|
{
|
|
idx_t i, j, nnbrs;
|
|
idx_t *nptr, *nind;
|
|
idx_t *xadj, *adjncy;
|
|
idx_t *marker, *nbrs;
|
|
|
|
|
|
/* construct the node-element list first */
|
|
nptr = ismalloc(nn+1, 0, "CreateGraphNodal: nptr");
|
|
nind = imalloc(eptr[ne], "CreateGraphNodal: nind");
|
|
|
|
for (i=0; i<ne; i++) {
|
|
for (j=eptr[i]; j<eptr[i+1]; j++)
|
|
nptr[eind[j]]++;
|
|
}
|
|
MAKECSR(i, nn, nptr);
|
|
|
|
for (i=0; i<ne; i++) {
|
|
for (j=eptr[i]; j<eptr[i+1]; j++)
|
|
nind[nptr[eind[j]]++] = i;
|
|
}
|
|
SHIFTCSR(i, nn, nptr);
|
|
|
|
|
|
/* Allocate memory for xadj, since you know its size.
|
|
These are done using standard malloc as they are returned
|
|
to the calling function */
|
|
if ((xadj = (idx_t *)malloc((nn+1)*sizeof(idx_t))) == NULL)
|
|
gk_errexit(SIGMEM, "***Failed to allocate memory for xadj.\n");
|
|
*r_xadj = xadj;
|
|
iset(nn+1, 0, xadj);
|
|
|
|
/* allocate memory for working arrays used by FindCommonElements */
|
|
marker = ismalloc(nn, 0, "CreateGraphNodal: marker");
|
|
nbrs = imalloc(nn, "CreateGraphNodal: nbrs");
|
|
|
|
for (i=0; i<nn; i++) {
|
|
xadj[i] = FindCommonNodes(i, nptr[i+1]-nptr[i], nind+nptr[i], eptr,
|
|
eind, marker, nbrs);
|
|
}
|
|
MAKECSR(i, nn, xadj);
|
|
|
|
/* Allocate memory for adjncy, since you now know its size.
|
|
These are done using standard malloc as they are returned
|
|
to the calling function */
|
|
if ((adjncy = (idx_t *)malloc(xadj[nn]*sizeof(idx_t))) == NULL) {
|
|
free(xadj);
|
|
*r_xadj = NULL;
|
|
gk_errexit(SIGMEM, "***Failed to allocate memory for adjncy.\n");
|
|
}
|
|
*r_adjncy = adjncy;
|
|
|
|
for (i=0; i<nn; i++) {
|
|
nnbrs = FindCommonNodes(i, nptr[i+1]-nptr[i], nind+nptr[i], eptr,
|
|
eind, marker, nbrs);
|
|
for (j=0; j<nnbrs; j++)
|
|
adjncy[xadj[i]++] = nbrs[j];
|
|
}
|
|
SHIFTCSR(i, nn, xadj);
|
|
|
|
gk_free((void **)&nptr, &nind, &marker, &nbrs, LTERM);
|
|
}
|
|
|
|
|
|
/*****************************************************************************/
|
|
/*! This function finds the union of nodes that are in the same elements with
|
|
the ``query'' node.
|
|
*/
|
|
/*****************************************************************************/
|
|
idx_t FindCommonNodes(idx_t qid, idx_t nelmnts, idx_t *elmntids, idx_t *eptr,
|
|
idx_t *eind, idx_t *marker, idx_t *nbrs)
|
|
{
|
|
idx_t i, ii, j, jj, k;
|
|
|
|
/* find all nodes that share at least one element with qid */
|
|
marker[qid] = 1; /* this is to prevent self-loops */
|
|
for (k=0, i=0; i<nelmnts; i++) {
|
|
j = elmntids[i];
|
|
for (ii=eptr[j]; ii<eptr[j+1]; ii++) {
|
|
jj = eind[ii];
|
|
if (marker[jj] == 0) {
|
|
nbrs[k++] = jj;
|
|
marker[jj] = 1;
|
|
}
|
|
}
|
|
}
|
|
|
|
/* reset the marker */
|
|
marker[qid] = 0;
|
|
for (i=0; i<k; i++) {
|
|
marker[nbrs[i]] = 0;
|
|
}
|
|
|
|
return k;
|
|
}
|
|
|
|
|
|
|
|
/*************************************************************************/
|
|
/*! This function creates and initializes a mesh_t structure */
|
|
/*************************************************************************/
|
|
mesh_t *CreateMesh(void)
|
|
{
|
|
mesh_t *mesh;
|
|
|
|
mesh = (mesh_t *)gk_malloc(sizeof(mesh_t), "CreateMesh: mesh");
|
|
|
|
InitMesh(mesh);
|
|
|
|
return mesh;
|
|
}
|
|
|
|
|
|
/*************************************************************************/
|
|
/*! This function initializes a mesh_t data structure */
|
|
/*************************************************************************/
|
|
void InitMesh(mesh_t *mesh)
|
|
{
|
|
memset((void *)mesh, 0, sizeof(mesh_t));
|
|
}
|
|
|
|
|
|
/*************************************************************************/
|
|
/*! This function deallocates any memory stored in a mesh */
|
|
/*************************************************************************/
|
|
void FreeMesh(mesh_t **r_mesh)
|
|
{
|
|
mesh_t *mesh = *r_mesh;
|
|
|
|
gk_free((void **)&mesh->eptr, &mesh->eind, &mesh->ewgt, &mesh, LTERM);
|
|
|
|
*r_mesh = NULL;
|
|
}
|
|
|
|
|