GALAHAD NODEND package#

purpose#

The nodend package finds a symmetric row and column permutation, $P A P^T$, of a symmetric, sparse matrix $A$ with the aim of limiting the fill-in during subsequent Cholesky-like factorization. The package is actually a wrapper to the METIS_NodeND procedure from versions 4.0, 5.1 and 5.2 of the METIS package from the Karypis Lab.

See Section 4 of $GALAHAD/doc/nodend.pdf for additional details.

method#

Variants of node-based nested-dissection ordering are employed.

reference#

The basic methods are those described in

G. Karypis. METIS, A software package for partitioning unstructured graphs, partitioning meshes, and computing fill-reducing orderings of sparse matrices, Version 5, Department of Computer Science & Engineering, University of Minnesota Minneapolis, MN 55455, USA (2013), see

KarypisLab/METIS

matrix storage#

The sparsity pattern of the symmetric $n$ by $n$ matrix $A$ may be presented and stored in a variety of formats. But crucially symmetry is exploited by only storing values from the lower triangular part (i.e, those entries that lie on or below the leading diagonal).

Sparse co-ordinate storage format: Only the nonzero entries of the matrices are stored. For the $l$-th entry, $1 \leq l \leq ne$, of $A$, its row index i and column index j, $1 \leq j \leq i \leq n$, are stored as the $l$-th components of the integer arrays A_row and A_col, respectively, while the number of nonzeros is recorded as A_ne = $ne$. Note that only the entries in the lower triangle should be stored. The string A_type = ‘coordinate’ should be specified.

Sparse row-wise storage format: Again only the nonzero entries are stored, but this time they are ordered so that those in row i appear directly before those in row i+1. For the i-th row of $A$ the i-th component of the integer array A_ptr holds the position of the first entry in this row, while A_ptr(n+1) holds the total number of entries plus one. The column indices j, $1 \leq j \leq i$, of the entries in the i-th row are stored in components l = A_ptr(i), …, A_ptr(i+1)-1 of the integer array A_col. Note that as before only the entries in the lower triangle should be stored. For sparse matrices, this scheme almost always requires less storage than its predecessor. The string A_type = ‘sparse_by_rows’ should be specified.

parametric integer type INT#

Below, the symbol INT refers to a parametric integer type that may be Int32 (32-bit integer) or Int64 (64-bit integer).

callable functions#

    function nodend_initialize(INT, data, control, status)

Set default control values and initialize private data

Parameters:

data

holds private internal data

control

is a structure containing control information (see nodend_control_type)

status

is a scalar variable of type INT that gives the exit status from the package. Possible values are (currently):

0
The initialization was successful.

    function nodend_read_specfile(INT, control, specfile)

Read the content of a specification file, and assign values associated with given keywords to the corresponding control parameters. An in-depth discussion of specification files is available, and a detailed list of keywords with associated default values is provided in $GALAHAD/src/nodend/NODEND.template. See also Table 2.1 in the Fortran documentation provided in $GALAHAD/doc/nodend.pdf for a list of how these keywords relate to the components of the control structure.

Parameters:

control	is a structure containing control information (see nodend_control_type)
specfile	is a one-dimensional array of type Vararg{Cchar} that must give the name of the specification file

    function nodend_order(INT, control, data, status, n, perm,
                          A_type, ne, A_row, A_col, A_ptr)

Order problem data into internal storage prior to solution.

Parameters:

control	is a structure whose members provide control parameters for the remaining procedures (see nodend_control_type)
data	holds private internal data
status	is a scalar variable of type INT that gives the exit status from the package. Possible values are: 1 The order was successful, and the package is ready for the solve phase -1 An allocation error occurred. A message indicating the offending array is written on unit control.error, and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively. -2 A deallocation error occurred. A message indicating the offending array is written on unit control.error and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively. -3 One of the restrictions n $> 0$, A.n $> 0$ or A.ne $< 0$, for co-ordinate entry, or requirements that A.type contain its relevant string ‘coordinate’ or ‘sparse_by_rows, and control.version in one of ‘4.0’, ‘5.1’ or ‘5.2’ has been violated. -26 The requested version of METIS is not available. -57 METIS has insufficient memory to continue. -71 An internal METIS error occurred.
n	is a scalar variable of type INT that holds the number of variables.
perm	is a one-dimensional array of size n and type INT, that returns the computed permutation array, so that the perm[i]-th rows and columns in the permuted matrix $P A P^T$ correspond to those labelled i in $A$, 0 $\leq$ i $\leq$ n-1.
A_type	is a one-dimensional array of type Vararg{Cchar} that specifies the symmetric storage scheme used $A$. It should be one of ‘coordinate’ or ‘sparse_by_rows’; lower or upper case variants are allowed. If A_type is not one of the supported values, the identity permutation will be returned.
ne	is a scalar variable of type INT that holds the number of entries in the lower triangular part of A in the sparse co-ordinate storage scheme. It need not be set for the other scheme.
A_row	is a one-dimensional array of size ne and type INT that holds the row indices of the lower triangular part of $A$ in the sparse co-ordinate storage scheme. It need not be set for the other scheme, and in this case can be C_NULL
A_col	is a one-dimensional array of size ne and type INT that holds the column indices of the lower triangular part of $A$.
A_ptr	is a one-dimensional array of size n+1 and type INT that holds the starting position of each row of the lower triangular part of $A$, as well as the total number of entries, in the sparse row-wise storage scheme. It need not be set when the other scheme is used, and in this case can be C_NULL

    function nodend_information(INT, data, inform, status)

Provides output information

Parameters:

data

holds private internal data

inform

is a structure containing output information (see nodend_inform_type)

status

is a scalar variable of type INT that gives the exit status from the package. Possible values are (currently):

0
The values were recorded successfully

available structures#

nodend_control_type structure#

    struct nodend_control_type{INT}
      f_indexing::Bool
      version::NTuple{31,Cchar}
      error::INT
      out::INT
      print_level::INT
      no_metis_4_use_5_instead::Bool
      prefix::NTuple{31,Cchar}
      metis4_ptype::INT
      metis4_ctype::INT
      metis4_itype::INT
      metis4_rtype::INT
      metis4_dbglvl::INT
      metis4_oflags::INT
      metis4_pfactor::INT
      metis4_nseps::INT
      metis5_ptype::INT
      metis5_objtype::INT
      metis5_ctype::INT
      metis5_iptype::INT
      metis5_rtype::INT
      metis5_dbglvl::INT
      metis5_niter::INT
      metis5_ncuts::INT
      metis5_seed::INT
      metis5_no2hop::INT
      metis5_minconn::INT
      metis5_contig::INT
      metis5_compress::INT
      metis5_ccorder::INT
      metis5_pfactor::INT
      metis5_nseps::INT
      metis5_ufactor::INT
      metis5_niparts::INT
      metis5_ondisk::INT
      metis5_dropedges::INT
      metis5_twohop::INT
      metis5_fast::INT

detailed documentation#

control derived type as a Julia structure

components#

Bool f_indexing

use C or Fortran sparse matrix indexing

char version[31]

specify the version of METIS to be used. Possible values are 4.0, 5.1 and 5.2

INT error

error and warning diagnostics occur on stream error

INT out

general output occurs on stream out

INT print_level

Bool no_metis_4_use_5_instead

if .no_metis_4_use_5_instead is true, and METIS 4.0 is not availble, use Metis 5.2 instead

NTuple{31,Cchar} prefix

all output lines will be prefixed by .prefix(2:LEN(TRIM(.prefix))-1) where .prefix contains the required string enclosed in quotes, e.g. “string” or ‘string’

INT metis4_ptype

the partitioning method employed. 0 = multilevel recursive bisectioning: 1 = multilevel k-way partitioning

INT metis4_ctype

the matching scheme to be used during coarsening: 1 = random matching, 2 = heavy-edge matching, 3 = sorted heavy-edge matching, and 4 = k-way sorted heavy-edge matching.

INT metis4_itype

the algorithm used during initial partitioning: 1 = edge-based region growing and 2 = node-based region growing.

INT metis4_rtype

the algorithm used for refinement: 1 = two-sided node Fiduccia-Mattheyses (FM) refinement, and 2 = one-sided node FM refinement.

INT metis4_dbglvl

the amount of progress/debugging information printed: 0 = nothing, 1 = timings, and $>$ 1 increasingly more.

INT metis4_oflags

select whether or not to compress the graph, and to order connected components separately: 0 = do neither, 1 = try to compress the graph, 2 = order each connected component separately, and 3 = do both.

INT metis4_pfactor

the minimum degree of the vertices that will be ordered last. More specifically, any vertices with a degree greater than 0.1 metis4_pfactor times the average degree are removed from the graph, an ordering of the rest of the vertices is computed, and an overall ordering is computed by ordering the removed vertices at the end of the overall ordering. Any value smaller than 1 means that no vertices will be ordered last.

INT metis4_nseps

the number of different separators that the algorithm will compute at each level of nested dissection.

INT metis5_ptype

the partitioning method. The value 0 gives multilevel recursive bisectioning, while 1 corresponds to multilevel $k$-way partitioning.

INT metis5_objtype

the type of the objective. Currently the only and default value metis5_objtype = 2, specifies node-based nested dissection, and any invalid value will be replaced by this default.

INT metis5_ctype

the matching scheme to be used during coarsening: 0 = random matching, and 1 = sorted heavy-edge matching.

INT metis5_iptype

the algorithm used during initial partitioning: 2 = derive separators from edge cuts, and 3 = grow bisections using a greedy node-based strategy.

INT metis5_rtype

the algorithm used for refinement: 2 = Two-sided node FM refinement, and 3 = One-sided node FM refinement.

INT metis5_dbglvl

the amount of progress/debugging information printed: 0 = nothing, 1 = diagnostics, 2 = plus timings, and $>$ 2 plus more.

INT metis5_niparts

the number of initial partitions used by MeTiS 5.2.

INT metis5_niter

the number of iterations used by the refinement algorithm.

INT metis5_ncuts

the number of different partitionings that it will compute: -1 = not used.

INT metis5_seed

the seed for the random number generator.

INT metis5_ondisk

whether on-disk storage is used (0 = no, 1 = yes) by MeTiS 5.2.

INT metis5_minconn

specify that the partitioning routines should try to minimize the maximum degree of the subdomain graph: 0 = no, 1 = yes, and -1 = not used.

INT metis5_contig

specify that the partitioning routines should try to produce partitions that are contiguous: 0 = no, 1 = yes, and -1 = not used.

INT metis5_compress

specify that the graph should be compressed by combining together vertices that have identical adjacency lists: 0 = no, and 1 = yes.

INT metis5_ccorder

specify if the connected components of the graph should first be identified and ordered separately: 0 = no, and 1 = yes.

INT metis5_pfactor

INT metis5_nseps

the number of different separators that the algorithm will compute at each level of nested dissection.

INT metis5_ufactor

the maximum allowed load imbalance (1 +metis5_ufactor)/1000 among the partitions.

INT metis5_dropedges

will edges be dropped (0 = no, 1 = yes) by MeTiS 5.2.

INT metis5_no2hop

specify that the coarsening will not perform any 2–hop matchings when the standard matching approach fails to sufficiently coarsen the graph: 0 = no, and 1 = yes.

INT metis5_twohop

reserved for future use but ignored at present.

INT metis5_fast

reserved for future use but ignored at present.

nodend_inform_type structure#

    struct nodend_inform_type{INT}
      status::INT
      alloc_status::INT
      bad_alloc::NTuple{81,Cchar}
      version::NTuple{4,Cchar}

detailed documentation#

inform derived type as a Julia structure

components#

INT status

return status. Possible values are:

0

The call was successful.
-1

An allocation error occurred. A message indicating the offending array is written on unit control.error, and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively.
-2

A deallocation error occurred. A message indicating the offending array is written on unit control.error and the returned allocation status and a string containing the name of the offending array are held in inform.alloc_status and inform.bad_alloc respectively.
-3

One of the restrictions n $> 0$, A.n $> 0$ or A.ne $< 0$, for co-ordinate entry, or requirements that A.type contain its relevant string ‘COORDINATE’ or ‘SPARSE_BY_ROWS’, and control.version in one of ‘4.0’, ‘5.1’ or ‘5.2’ has been violated.
-26

The requested version of METIS is not available.
-57

METIS has insufficient memory to continue.
-71

An internal METIS error occurred.

INT alloc_status

the status of the last attempted allocation/deallocation

NTuple{81,Cchar} bad_alloc

the name of the array for which an allocation/deallocation error occurred

NTuple{4,Cchar} version

the version of METIS that was actually used

Index

example calls#

This is an example of how to use the package to pick row/column orderings of a given matrix $A$ prior to Cholesky-like factorization; the code is available in $GALAHAD/src/nodend/Julia/test_nodend.jl . A variety of supported matrix storage formats are shown.

# test_nodend.jl
# Simple code to test the Julia interface to NODEND

using GALAHAD
using Test
using Printf
using Accessors
using Quadmath

function test_nodend(::Type{INT}) where {INT}
  # Derived types
  data = Ref{Ptr{Cvoid}}()
  control = Ref{nodend_control_type{INT}}()
  inform = Ref{nodend_inform_type{INT}}()

  # Set problem data
  n = INT(10)  # dimension
  A_ne = 2 * n - 1 # Hesssian elements, NB lower triangle
  A_dense_ne = div(n * (n + 1), 2) # dense Hessian elements
  A_row = zeros(INT, A_ne) # row indices,
  A_col = zeros(INT, A_ne) # column indices
  A_ptr = zeros(INT, n + 1)  # row pointers

  # Set output storage
  perm = zeros(INT, n) # permutation vector
  st = ' '
  status = Ref{INT}()

  # A = tridiag(2,1)
  l = 1
  A_ptr[1] = l
  A_row[l] = 1
  A_col[l] = 1
  for i in 2:n
    l = l + 1
    A_ptr[i] = l
    A_row[l] = i
    A_col[l] = i - 1
    l = l + 1
    A_row[l] = i
    A_col[l] = i
  end
  A_ptr[n + 1] = l + 1

  @printf(" fortran sparse matrix indexing\n\n")
  @printf(" basic tests of nodend storage formats\n\n")

  for d in 1:3

    # Initialize NODEND
    nodend_initialize(INT, data, control, status)

    # Set user-defined control options
    @reset control[].f_indexing = true # fortran sparse matrix indexing

    # sparse co-ordinate storage
    if d == 1
      st = 'C'
      nodend_order(INT, control, data, status, n, perm,
                   "coordinate", A_ne, A_row, A_col, C_NULL)
    end

    # sparse by rows
    if d == 2
      st = 'R'
      nodend_order(INT, control, data, status, n, perm,
                   "sparse_by_rows", A_ne, C_NULL, A_col, A_ptr)
    end

    # dense
    if d == 3
      st = 'D'
      nodend_order(INT, control, data, status, n, perm,
                   "dense", A_dense_ne, C_NULL, C_NULL, C_NULL)
    end

    nodend_information(INT, data, inform, status)

    if inform[].status == 0
      @printf("%c: NODEND_order success, perm: ", st);
      for i in 1:n
        @printf("%1i ", perm[i]);
      end
      @printf("\n");
    else
      @printf("%c: NODEND_solve exit status = %1i\n", st, inform[].status)
    end

    # Delete internal workspace
    nodend_terminate(INT, data, control, inform)
  end
end

for (INT, libgalahad) in ((Int32, GALAHAD.libgalahad_single      ),
                          (Int64, GALAHAD.libgalahad_single_64   ),
                          (Int32, GALAHAD.libgalahad_double      ),
                          (Int64, GALAHAD.libgalahad_double_64   ),
                          (Int32, GALAHAD.libgalahad_quadruple   ),
                          (Int64, GALAHAD.libgalahad_quadruple_64))
  if isfile(libgalahad)
    @testset "NODEND -- $INT" begin
      @test test_nodend(INT) == 0
    end
  end
end