METIS_NodeND Interface

interface
public function METIS_NodeND(nvtxs, xadj, adjncy, vwgt, options, perm, iperm) result(ierr) bind(C,name="METIS_NodeND")

Arguments

Type IntentOptional AttributesName
integer(kind=idx_t), intent(in) :: nvtxs

The number of vertices in the graph.

integer(kind=idx_t), intent(in), dimension(*):: xadj

The adjacency structure of the graph as described in Section 5.5 of the manual.

integer(kind=idx_t), intent(in), dimension(*):: adjncy

The adjacency structure of the graph as described in Section 5.5 of the manual.

integer(kind=idx_t), intent(in), optional :: vwgt(nvtxs)

An array of size nvtxs specifying the weights of the vertices.

integer(kind=idx_t), intent(in), optional :: options(METIS_NOPTIONS)

This is the array of options as described in Section 5.4 of the manual. See description for valid options.

integer(kind=idx_t), intent(out) :: perm(nvtxs)

Vectors of size nvtxs. Upon successful completion, they store the fill-reducing permutation and inverse-permutation. More in the description.

integer(kind=idx_t), intent(out) :: iperm(nvtxs)

Vectors of size nvtxs. Upon successful completion, they store the fill-reducing permutation and inverse-permutation. More in the description.

Return Value integer(kind=idx_t)

METIS_OK - Indicates that the function returned normally.
METIS_ERROR_INPUT - Indicates an input error.
METIS_ERROR_MEMORY - Indicates that it could not allocate the required memory.
METIS_ERROR - Indicates some other type of error.

Description

This function computes fill reducing orderings of sparse matrices using the multilevel nested dissection algorithm.

Let \(A\) be the original matrix and \(A'\) be the permuted matrix. The arrays perm and iperm are defined as follows. Row (column) i of \(A'\) is the perm(i) row (column) of \(A\), and row (column) i of \(A\) is the iperm(i) row (column) of \(A'\). the numbering of this vector starts from either 0 or 1, depending on the value of options(METIS_OPTION_NUMBERING).

If the graph is weighted, meaning vgwt was provided, the nested dissection ordering computes vertex separators that minimize the sum of the weights of the vertices on the separators.

The following options are valid:
METIS_OPTION_CTYPE, METIS_OPTION_RTYPE, METIS_OPTION_NO2HOP, METIS_OPTION_NSEPS, METIS_OPTION_NITER, METIS_OPTION_UFACTOR, METIS_OPTION_COMPRESS, METIS_OPTION_CCORDER, METIS_OPTION_SEED, METIS_OPTION_PFACTOR, METIS_OPTION_NUMBERING, METIS_OPTION_DBGLVL

Example

The code below generates the nodal graph of the following mesh:

 use metis_interface, only: idx_t, METIS_SetDefaultOptions, METIS_NodeND
integer(kind=idx_t), parameter :: n = 15, m = 22
integer(kind=idx_t) :: xadj(n+1), adjncy(2*m)
integer(kind=idx_t) :: perm(n), iperm(n)
integer(kind=idx_t) :: options(0:39), ierr

xadj = [1,3,6,9,12,14,17,21,25,29,32,34,37,40,43,45]
adjncy = [2,6,1,3,7,2,4,8,3,5,9,4,10,1,7,11,2,6, &
                  8,12,3,7,9,13,4,8,10,14,5,9,15,6,12,7,11,13, &
                  8,12,14,9,13,15,10,14]

ierr = METIS_SetDefaultOptions(options)
options(18) = 1

ierr = METIS_NodeND(n,xadj,adjncy,options=options,perm=perm,iperm=iperm)
end