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$ sparse matrix $A$ may be presented and stored in a couple of formats. But crucially symmetry is exploited by only indices 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, $0 \leq l \leq ne-1$, of $H$, its row index i and column index j, $0 \leq j \leq i \leq n-1$, are stored as the $l$-th components of the integer arrays H_row and H_col, respectively, while the number of nonzeros is recorded as H_ne = $ne$. Note that only the entries in the lower triangle should be stored. The string H_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 $H$ the i-th component of the integer array H_ptr holds the position of the first entry in this row, while H_ptr(n) holds the total number of entries. The column indices j, $0 \leq j \leq i$, are stored in components l = H_ptr(i), …, H_ptr(i+1)-1 of the integer array H_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 H_type = ‘sparse_by_rows’ should be specified.

introduction to function calls#

To find the required permutation, functions from the nodend package must be called in the following order:

nodend_initialize - provide default control parameters and set up initial data structures
nodend_read_specfile (optional) - override control values by reading replacement values from a file
nodend_import find a row/colum permutation for sparse Cholesky-like factorization
nodend_information (optional) - recover information about the solution and solution process

See the examples section for illustrations of use.

callable functions#

overview of functions provided#

// namespaces

namespace conf;

// typedefs

typedef double rpc_;
typedef int ipc_;

// structs

struct nodend_control_type;
struct nodend_time_type;
struct nodend_inform_type;

// global functions

void nodend_initialize(
    void **data,
    struct nodend_control_type* control,
    ipc_ *status
);

void nodend_read_specfile(struct nodend_control_type* control, const char specfile[]);

void nodend_order(
    struct nodend_control_type* control,
    void **data,
    ipc_ *status,
    ipc_ n,
    ipc_ perm[],
    const char H_type[],
    ipc_ ne,
    const ipc_ A_row[],
    const ipc_ A_col[],
    const ipc_ A_ptr[]
);

void nodend_information(
    void **data,
    struct nodend_inform_type* inform,
    ipc_ *status
);

void nodend_terminate(
    void **data
);

typedefs#

typedef double rpc_

rpc_ is the real working precision used, but may be changed to float by defining the preprocessor variable REAL_32 or (if supported) to __real128 using the variable REAL_128.

typedef int ipc_

ipc_ is the default integer word length used, but may be changed to int64_t by defining the preprocessor variable INTEGER_64.

function calls#

void nodend_initialize(
    void **data,
    struct nodend_control_type* control,
    ipc_ *status
)

Set default control values and initialize private data

Parameters:

data

holds private internal data

control

is a struct containing control information (see nodend_control_type)

status

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

0
The initialization was successful.

void nodend_read_specfile(struct nodend_control_type* control, const char 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 struct containing control information (see nodend_control_type)
specfile	is a character string containing the name of the specification file

void nodend_order(
    struct nodend_control_type* control,
    void **data,
    ipc_ *status,
    ipc_ n,
    ipc_ perm[],
    const char H_type[],
    ipc_ ne,
    const ipc_ A_row[],
    const ipc_ A_col[],
    const ipc_ A_ptr[]
)

Find a row/colum permutation for sparse Cholesky-like factorization.

Parameters:

control	is a struct whose members provide control paramters for the remaining prcedures (see nodend_control_type)
data	holds private internal data
status	is a scalar variable of type ipc_, that gives the exit status from the package. Possible values are: 1 The import 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 The restriction n > 0 or requirement that type contains its relevant string ‘coordinate’ or ‘sparse_by_rows’ has been violated. 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 ipc_, that holds the number of variables.
perm	is a one-dimensional array of size n and type ipc_, 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 char that specifies the symmetric storage scheme used for the Hessian. 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 ipc_, that holds the number of entries in the lower triangular part of H in the sparse co-ordinate storage scheme. It need not be set for any of the other three schemes.
A_row	is a one-dimensional array of size ne and type ipc_, 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 NULL
A_col	is a one-dimensional array of size ne and type ipc_, 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 ipc_, 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 NULL

void nodend_information(
    void **data,
    struct nodend_inform_type* inform,
    ipc_ *status
)

Provides output information

Parameters:

data

holds private internal data

inform

is a struct containing output information (see nodend_inform_type)

status

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

0
The values were recorded successfully

void nodend_terminate(
    void **data
)

Deallocate all internal private storage

Parameters:

data	holds private internal data

available structures#

nodend_control_type structure#

#include <galahad_nodend.h>

struct nodend_control_type {
  // fields

  bool f_indexing;
  char version[31];
  ipc_ error;
  ipc_ out;
  ipc_ print_level;
  bool no_metis_4_use_5_instead;
  char prefix[31];
  ipc_ metis4_ptype;
  ipc_ metis4_ctype;
  ipc_ metis4_itype;
  ipc_ metis4_rtype;
  ipc_ metis4_dbglvl;
  ipc_ metis4_oflags;
  ipc_ metis4_pfactor;
  ipc_ metis4_nseps;
  ipc_ metis5_ptype;
  ipc_ metis5_objtype;
  ipc_ metis5_ctype;
  ipc_ metis5_iptype;
  ipc_ metis5_rtype;
  ipc_ metis5_dbglvl;
  ipc_ metis5_niter;
  ipc_ metis5_ncuts;
  ipc_ metis5_seed;
  ipc_ metis5_no2hop;
  ipc_ metis5_minconn;
  ipc_ metis5_contig;
  ipc_ metis5_compress;
  ipc_ metis5_ccorder;
  ipc_ metis5_pfactor;
  ipc_ metis5_nseps;
  ipc_ metis5_ufactor;
  ipc_ metis5_niparts;
  ipc_ metis5_ondisk;
  ipc_ metis5_dropedges;
  ipc_ metis5_twohop;
  ipc_ metis5_fast;
  };

detailed documentation#

control derived type as a C struct

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

ipc_ error

error and warning diagnostics occur on stream error

ipc_ out

general output occurs on stream out

ipc_ print_level

the level of output required is specified by 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

char prefix[31]

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’

ipc_ metis4_ptype

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

ipc_ 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.

ipc_ metis4_itype

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

ipc_ metis4_rtype

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

ipc_ metis4_dbglvl

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

ipc_ 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.

ipc_ 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.

ipc_ metis4_nseps

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

ipc_ metis5_ptype

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

ipc_ 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.

ipc_ metis5_ctype

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

ipc_ metis5_iptype

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

ipc_ metis5_rtype

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

ipc_ metis5_dbglvl

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

ipc_ metis5_niparts

the number of initial partitions used by MeTiS 5.2.

ipc_ metis5_niter

the number of iterations used by the refinement algorithm.

ipc_ metis5_ncuts

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

ipc_ metis5_seed

the seed for the random number generator.

ipc_ metis5_ondisk

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

ipc_ 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.

ipc_ metis5_contig

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

ipc_ metis5_compress

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

ipc_ metis5_ccorder

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

ipc_ metis5_pfactor

ipc_ metis5_nseps

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

ipc_ metis5_ufactor

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

ipc_ metis5_dropedges

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

ipc_ 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.

ipc_ metis5_twohop

reserved for future use but ignored at present.

ipc_ metis5_fast

reserved for future use but ignored at present.

nodend_time_type structure#

#include <galahad_nodend.h>

struct nodend_time_type {
    // components

    rpc_ total;
    rpc_ metis;
    rpc_ clock_total;
    rpc_ clock_metis;
};

detailed documentation#

time derived type as a C struct

components#

rpc_ total

the total CPU time spent in the package

rpc_ metis

the CPU time spent in the METIS package

rpc_ clock_total

the total clock time spent in the package

rpc_ clock_metis

the clock time spent in the METIS package

nodend_inform_type structure#

#include <galahad_nodend.h>

struct nodend_inform_type {
  // fields

  ipc_ status;
  ipc_ alloc_status;
  char bad_alloc[81];
  char version[4];
  struct nodend_time_type time;
  };

detailed documentation#

inform derived type as a C struct

components#

ipc_ 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’, ‘SPARSE_BY_ROWS’ or ‘DENSE’, 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.

ipc_ alloc_status

the status of the last attempted allocation/deallocation

char bad_alloc[81]

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

char version[4]

the version of METIS that was actually used

struct nodend_time_type time

timings (see above)

Index

example calls#

This is an example of how to use the package to find a suitable row/column permutation of a given matrix; the code is available in $GALAHAD/src/nodend/C/nodendt.c .

Notice that C-style indexing is used, and that this is flagged by setting control.f_indexing to false. The integer type ipc_ from galahad_precision.h is set to int by default, but to int64_t if the preprocessor variable INTEGER_64 is defined.

/* nodendt.c */
/* Full test for the NODEND C interface using C sparse matrix indexing */

#include <stdio.h>
#include <math.h>
#include "galahad_precision.h"
#include "galahad_cfunctions.h"
#include "galahad_nodend.h"

int main(void) {

    // Derived types
    void *data;
    struct nodend_control_type control;
    struct nodend_inform_type inform;

    // Set problem data
    ipc_ n = 10; // dimension
    ipc_ A_ne = 2 * n - 1; // Hesssian elements, NB lower triangle
    ipc_ A_dense_ne = n * ( n + 1 ) / 2; // dense Hessian elements
    ipc_ A_row[A_ne]; // row indices,
    ipc_ A_col[A_ne]; // column indices
    ipc_ A_ptr[n+1];  // row pointers
    ipc_ perm[n]; // permutation,

    // Set output storage
    char st = ' ';
    ipc_ i, l, status;

    // A = tridiag(2,1)

    l = 0 ;
    A_ptr[0] = l;
    A_row[l] = 0; A_col[l] = 0;
    for( ipc_ i = 1; i < n; i++)
    {
      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;
    }
    A_ptr[n] = l + 1;

    printf(" C sparse matrix indexing\n\n");

    printf(" basic tests of nodend storage formats\n\n");

    for( ipc_ d=1; d <= 3; d++){

        // Initialize NODEND
        nodend_initialize( &data, &control, &status );

        // Set user-defined control options
        control.f_indexing = false; // C sparse matrix indexing

        switch(d){
            case 1: // sparse co-ordinate storage
                st = 'C';
                nodend_order( &control, &data, &status, n, perm,
                              "coordinate", A_ne, A_row, A_col, NULL );
                break;
            printf(" case %1" i_ipc_ " break\n",d);
            case 2: // sparse by rows
                st = 'R';
                nodend_order( &control, &data, &status, n, perm,
                              "sparse_by_rows", A_ne, NULL, A_col, A_ptr );
                break;
            case 3: // dense
                st = 'D';
                nodend_order( &control, &data, &status, n, perm,
                              "dense", A_dense_ne, NULL, NULL, NULL );
                break;
            }
        nodend_information( &data, &inform, &status );

        if(inform.status == 0){
            printf("%c: NODEND_order success, perm: ", st);
            for( i = 0; i < n; i++) printf("%1" i_ipc_ " ", perm[i]);
            printf("\n");
        }else{
            printf("%c: NODEND_order exit status = %1" i_ipc_ "\n",
                   st, inform.status);
        }

        // Terminate NODEND
        nodend_terminate( &data );

    }
}

This is the same example, but now fortran-style indexing is used; the code is available in $GALAHAD/src/nodend/C/nodendtf.c .

/* nodendtf.c */
/* Full test for the NODEND C interface using fortran sparse matrix indexing */

#include <stdio.h>
#include <math.h>
#include "galahad_precision.h"
#include "galahad_cfunctions.h"
#include "galahad_nodend.h"

int main(void) {

    // Derived types
    void *data;
    struct nodend_control_type control;
    struct nodend_inform_type inform;

    // Set problem data
    ipc_ n = 10; // dimension
    ipc_ A_ne = 2 * n - 1; // Hesssian elements, NB lower triangle
    ipc_ A_dense_ne = n * ( n + 1 ) / 2; // dense Hessian elements
    ipc_ A_row[A_ne]; // row indices,
    ipc_ A_col[A_ne]; // column indices
    ipc_ A_ptr[n+1];  // row pointers
    ipc_ perm[n]; // permutation,

    // Set output storage
    char st = ' ';
    ipc_ i, l, status;

    // A = tridiag(2,1)

    l = 0 ;
    A_ptr[0] = l + 1;
    A_row[l] = 1; A_col[l] = 1;
    for( ipc_ i = 1; i < n; i++)
    {
      l = l + 1;
      A_ptr[i] = l + 1;
      A_row[l] = i + 1; A_col[l] = i;
      l = l + 1;
      A_row[l] = i + 1; A_col[l] = i + 1;
    }
    A_ptr[n] = l + 2;

    printf(" Fortran sparse matrix indexing\n\n");

    printf(" basic tests of nodend storage formats\n\n");

    for( ipc_ d=1; d <= 3; d++){

        // Initialize NODEND
        nodend_initialize( &data, &control, &status );

        // Set user-defined control options
        control.f_indexing = true; // Fortran sparse matrix indexing

        switch(d){
            case 1: // sparse co-ordinate storage
                st = 'C';
                nodend_order( &control, &data, &status, n, perm,
                              "coordinate", A_ne, A_row, A_col, NULL );
                break;
            printf(" case %1" i_ipc_ " break\n",d);
            case 2: // sparse by rows
                st = 'R';
                nodend_order( &control, &data, &status, n, perm,
                              "sparse_by_rows", A_ne, NULL, A_col, A_ptr );
                break;
            case 3: // dense
                st = 'D';
                nodend_order( &control, &data, &status, n, perm,
                              "dense", A_dense_ne, NULL, NULL, NULL );
                break;
            }
        nodend_information( &data, &inform, &status );

        if(inform.status == 0){
            printf("%c: NODEND_order success, perm: ", st);
            for( i = 0; i < n; i++) printf("%1" i_ipc_ " ", perm[i]);
            printf("\n");
        }else{
            printf("%c: NODEND_order exit status = %1" i_ipc_ "\n",
                   st, inform.status);
        }

        // Terminate NODEND
        nodend_terminate( &data );

    }
}