The namespace for ALP/GraphBLAS algorithms. More...

Namespaces
	pregel
	The namespace for ALP/Pregel algorithms.

Classes
class	matrices
	Factories for creating ALP/GraphBLAS matrices with standard patterns such as identity matrices. More...

class	matrices< void, mode, backend, RIT, CIT, NIT >
	Factories for creating ALP/GraphBLAS pattern matrices. More...

Functions
template<Descriptor descr = descriptors::no_operation, typename IOType , typename NonzeroType , typename InputType , typename ResidualType , class Semiring = Semiring< operators::add< InputType, InputType, InputType >, operators::mul< IOType, NonzeroType, InputType >, identities::zero, identities::one >, class Minus = operators::subtract< ResidualType >, class Divide = operators::divide< ResidualType >>
RC	bicgstab (grb::Vector< IOType > &x, const grb::Matrix< NonzeroType > &A, const grb::Vector< InputType > &b, const size_t max_iterations, ResidualType tol, size_t &iterations, ResidualType &residual, Vector< InputType > &r, Vector< InputType > &rhat, Vector< InputType > &p, Vector< InputType > &v, Vector< InputType > &s, Vector< InputType > &t, const Semiring &semiring=Semiring(), const Minus &minus=Minus(), const Divide &divide=Divide())
	Solves a linear system \( b = Ax \) with \( x \) unknown by using the bi-conjugate gradient (bi-CG) stabilised method; i.e., BiCGstab. More...

template<Descriptor descr = descriptors::no_operation, typename IOType , typename ResidualType , typename NonzeroType , typename InputType , class Ring = Semiring< grb::operators::add< IOType >, grb::operators::mul< IOType >, grb::identities::zero, grb::identities::one >, class Minus = operators::subtract< IOType >, class Divide = operators::divide< IOType >, typename RSI , typename NZI , Backend backend>
grb::RC	conjugate_gradient (grb::Vector< IOType, backend > &x, const grb::Matrix< NonzeroType, backend, RSI, RSI, NZI > &A, const grb::Vector< InputType, backend > &b, const size_t max_iterations, ResidualType tol, size_t &iterations, ResidualType &residual, grb::Vector< IOType, backend > &r, grb::Vector< IOType, backend > &u, grb::Vector< IOType, backend > &temp, const Ring &ring=Ring(), const Minus &minus=Minus(), const Divide &divide=Divide())
	Solves a linear system \( b = Ax \) with \( x \) unknown by the Conjugate Gradients (CG) method on general fields. More...

template<Descriptor descr = descriptors::no_operation, typename OutputType , typename InputType1 , typename InputType2 , class Ring , class Division = grb::operators::divide< typename Ring::D3, typename Ring::D3, typename Ring::D4 >>
RC	cosine_similarity (OutputType &similarity, const Vector< InputType1 > &x, const Vector< InputType2 > &y, const Ring &ring=Ring(), const Division &div=Division())
	Computes the cosine similarity. More...

template<Descriptor descr = descriptors::no_operation, typename NonzeroType , typename ResidualType , class Ring = Semiring< grb::operators::add< NonzeroType >, grb::operators::mul< NonzeroType >, grb::identities::zero, grb::identities::one >, class Minus = operators::subtract< NonzeroType >, class Divide = operators::divide< NonzeroType >>
grb::RC	gmres (grb::Vector< NonzeroType > &x, const grb::Matrix< NonzeroType > &A, const grb::Vector< NonzeroType > &b, const size_t n_restart, const size_t max_iterations, const ResidualType tol, size_t &iterations, size_t &iterations_gmres, size_t &iterations_arnoldi, ResidualType &residual, std::vector< grb::Vector< NonzeroType > > &Q, std::vector< NonzeroType > &Hmatrix, grb::Vector< NonzeroType > &temp, std::vector< NonzeroType > &temp3, const Ring &ring=Ring(), const Minus &minus=Minus(), const Divide &divide=Divide(), const std::function< ResidualType(ResidualType) > &sqrtX=std_sqrt< ResidualType, ResidualType >)
	Solves a linear system \( b = Ax \) with \( x \) unknown using the Generalised Minimal Residual (GMRES) method on general fields. More...

template<typename NonzeroType , typename DimensionType , typename ResidualType , class Ring = Semiring< grb::operators::add< NonzeroType >, grb::operators::mul< NonzeroType >, grb::identities::zero, grb::identities::one >, class Minus = operators::subtract< NonzeroType >, class Divide = operators::divide< NonzeroType >>
grb::RC	hessolve (std::vector< NonzeroType > &H, const DimensionType n, const DimensionType &kspspacesize, const ResidualType tol, std::vector< NonzeroType > &vecx, const Ring &ring=Ring(), const Minus &minus=Minus(), const Divide &divide=Divide(), const std::function< ResidualType(ResidualType) > &sqrtX=std_sqrt< ResidualType, ResidualType >)
	Solves the least linear square problem of size n - 1, defined by the equation \( A x = b, \) where the \( A \) is an upper-Hessenberg matrix of size n x kspspacesize, and the vector \( b \), of length n, has only the first element nonzero. More...

template<Descriptor descr = descriptors::no_operation, bool criticalSection = false, typename IOType , typename NZType >
RC	kcore_decomposition (const Matrix< NZType > &A, Vector< IOType > &core, Vector< IOType > &distances, Vector< IOType > &temp, Vector< IOType > &update, Vector< bool > &status, IOType &k)
	The \( k \)-core decomposition algorithm. More...

template<Descriptor descr = descriptors::no_operation, typename IOType = double, class Operator = operators::square_diff< IOType, IOType, IOType >>
RC	kmeans_iteration (Matrix< IOType > &K, Vector< std::pair< size_t, IOType > > &clusters_and_distances, const Matrix< IOType > &X, const size_t max_iter=1000, const Operator &dist_op=Operator())
	The kmeans iteration given an initialisation. More...

template<Descriptor descr, typename OutputType , typename InputType >
RC	knn (Vector< OutputType > &u, const Matrix< InputType > &A, const size_t source, const size_t k, Vector< bool > &buf1)
	Given a graph and a source vertex, indicates which vertices are contained within k hops. More...

template<Descriptor descr = descriptors::no_operation, typename IOType = double, class Operator = operators::square_diff< IOType, IOType, IOType >>
RC	kpp_initialisation (Matrix< IOType > &K, const Matrix< IOType > &X, const Operator &dist_op=Operator())
	a simple implementation of the k++ initialisation algorithm for kmeans More...

template<typename IOType >
RC	label (Vector< IOType > &out, const Vector< IOType > &y, const Matrix< IOType > &W, const size_t n, const size_t l, const size_t maxIterations=1000)
	The label propagation algorithm. More...

template<Descriptor descr, class Ring , typename IOType , typename InputType >
RC	mpv (Vector< IOType > &u, const Matrix< InputType > &A, const size_t k, const Vector< IOType > &v, Vector< IOType > &temp, const Ring &ring)
	The matrix powers kernel. More...

template<Descriptor descr = descriptors::no_operation, class Ring , typename InputType , typename OutputType , Backend backend, typename Coords >
RC	norm2 (OutputType &x, const Vector< InputType, backend, Coords > &y, const Ring &ring=Ring(), const std::function< OutputType(OutputType) > sqrtX=std_sqrt< OutputType, OutputType >, const typename std::enable_if< std::is_floating_point< OutputType >::value, void >::type *=nullptr)
	Provides a generic implementation of the 2-norm computation. More...

template<Descriptor descr = descriptors::no_operation, bool preconditioned = true, typename IOType , typename ResidualType , typename NonzeroType , typename InputType , class Ring = Semiring< grb::operators::add< IOType >, grb::operators::mul< IOType >, grb::identities::zero, grb::identities::one >, class Minus = operators::subtract< IOType >, class Divide = operators::divide< IOType >, typename RSI , typename NZI , Backend backend>
grb::RC	preconditioned_conjugate_gradient (grb::Vector< IOType, backend > &x, const grb::Matrix< NonzeroType, backend, RSI, RSI, NZI > &A, const grb::Vector< InputType, backend > &b, const std::function< grb::RC(grb::Vector< IOType, backend > &, const grb::Vector< IOType, backend > &) > &Minv, const size_t max_iterations, ResidualType tol, size_t &iterations, ResidualType &residual, grb::Vector< IOType, backend > &r, grb::Vector< IOType, backend > &u, grb::Vector< IOType, backend > &temp, grb::Vector< IOType, backend > &temp_precond, const Ring &ring=Ring(), const Minus &minus=Minus(), const Divide &divide=Divide())
	Solves a preconditioned linear system \( b = M^{-1}Ax \) with \( x \) unknown by the Conjugate Gradients (CG) method on general fields. More...

template<Descriptor descr = descriptors::no_operation, typename NonzeroType , typename ResidualType , class Ring = Semiring< grb::operators::add< NonzeroType >, grb::operators::mul< NonzeroType >, grb::identities::zero, grb::identities::one >, class Minus = operators::subtract< NonzeroType >, class Divide = operators::divide< NonzeroType >>
grb::RC	preconditioned_gmres (grb::Vector< NonzeroType > &x, const grb::Matrix< NonzeroType > &M, const grb::Matrix< NonzeroType > &A, const grb::Vector< NonzeroType > &b, const size_t n_restart, const size_t max_iterations, const ResidualType tol, size_t &iterations, size_t &iterations_gmres, size_t &iterations_arnoldi, ResidualType &residual, std::vector< grb::Vector< NonzeroType > > &Q, std::vector< NonzeroType > &Hmatrix, grb::Vector< NonzeroType > &temp, std::vector< NonzeroType > &temp3, const Ring &ring=Ring(), const Minus &minus=Minus(), const Divide &divide=Divide(), const std::function< ResidualType(ResidualType) > &sqrtX=std_sqrt< ResidualType, ResidualType >)
	Solves a linear system \( b = Ax \) with \( x \) unknown using the Generalised Minimal Residual (GMRES) method on general fields. More...

template<Descriptor descr = descriptors::no_operation, typename IOType , typename NonzeroT >
RC	simple_pagerank (Vector< IOType > &pr, const Matrix< NonzeroT > &L, Vector< IOType > &pr_next, Vector< IOType > &pr_nextnext, Vector< IOType > &row_sum, const IOType alpha=0.85, const IOType conv=0.0000001, const size_t max=1000, size_t const iterations=nullptr, double const quality=nullptr)
	The canonical PageRank algorithm. More...

template<Descriptor descr = descriptors::no_operation, typename IOType , typename WeightType , typename BiasType , class ReluMonoid = Monoid< grb::operators::relu< IOType >, grb::identities::negative_infinity >, class Ring = Semiring< grb::operators::add< IOType >, grb::operators::mul< IOType >, grb::identities::zero, grb::identities::one >>
grb::RC	sparse_nn_single_inference (grb::Vector< IOType > &out, const grb::Vector< IOType > &in, const std::vector< grb::Matrix< WeightType > > &layers, const std::vector< BiasType > &biases, grb::Vector< IOType > &temp, const ReluMonoid &relu=ReluMonoid(), const Ring &ring=Ring())
	Performs an inference step of a single data element through a Sparse Neural Network defined by num_layers sparse weight matrices and num_layers biases. More...

template<Descriptor descr = descriptors::no_operation, typename IOType , typename WeightType , typename BiasType , typename ThresholdType = IOType, class MinMonoid = Monoid< grb::operators::min< IOType >, grb::identities::infinity >, class ReluMonoid = Monoid< grb::operators::relu< IOType >, grb::identities::negative_infinity >, class Ring = Semiring< grb::operators::add< IOType >, grb::operators::mul< IOType >, grb::identities::zero, grb::identities::one >>
grb::RC	sparse_nn_single_inference (grb::Vector< IOType > &out, const grb::Vector< IOType > &in, const std::vector< grb::Matrix< WeightType > > &layers, const std::vector< BiasType > &biases, const ThresholdType threshold, grb::Vector< IOType > &temp, const ReluMonoid &relu=ReluMonoid(), const MinMonoid &min=MinMonoid(), const Ring &ring=Ring())
	Performs an inference step of a single data element through a Sparse Neural Network defined by num_layers sparse weight matrices and num_layers biases. More...

template<bool normalize = false, typename IOType , typename InputType >
RC	spy (grb::Matrix< IOType > &out, const grb::Matrix< InputType > &in)
	Given an input matrix and a smaller output matrix, map nonzeroes from the input matrix into the smaller one and count the number of nonzeroes that are mapped from the bigger matrix into the smaller. More...

template<bool normalize = false, typename IOType >
RC	spy (grb::Matrix< IOType > &out, const grb::Matrix< bool > &in)
	Specialisation for boolean input matrices in. More...

template<bool normalize = false, typename IOType >
RC	spy (grb::Matrix< IOType > &out, const grb::Matrix< void > &in)
	Specialisation for void input matrices in. More...

template<typename OutputType , typename InputType >
OutputType	std_sqrt (const InputType x)
	An alias of std::sqrt where the input and output types are templated separately. More...

Detailed Description

The namespace for ALP/GraphBLAS algorithms.

Function Documentation

◆ bicgstab()

RC grb::algorithms::bicgstab	(	grb::Vector< IOType > &	x,
		const grb::Matrix< NonzeroType > &	A,
		const grb::Vector< InputType > &	b,
		const size_t	max_iterations,
		ResidualType	tol,
		size_t &	iterations,
		ResidualType &	residual,
		Vector< InputType > &	r,
		Vector< InputType > &	rhat,
		Vector< InputType > &	p,
		Vector< InputType > &	v,
		Vector< InputType > &	s,
		Vector< InputType > &	t,
		const Semiring &	semiring = `Semiring()`,
		const Minus &	minus = `Minus()`,
		const Divide &	divide = `Divide()`
	)

Solves a linear system \( b = Ax \) with \( x \) unknown by using the bi-conjugate gradient (bi-CG) stabilised method; i.e., BiCGstab.

Template Parameters

descr	Any descriptor to use for the computation (optional).
IOType	The solution vector element type.
NonzeroType	The system matrix entry type.
InputType	The element type of the right-hand side vector.
ResidualType	The type of the residuals used during computation.
Semiring	The semiring under which to perform the BiCGstab
Minus	The inverse operator of the additive operator of Semiring.
Divide	The inverse of the multiplicative operator of Semiring.

By default, these will be the regular add, mul, subtract, and divide over the types IOType, NonzeroType, InputType, and/or ResidualType, as appropriate.

Does not perform any preconditioning.

Parameters

[in,out]	x	On input: an initial guess to the solution \( Ax = b \). On output: if grb::SUCCESS is returned, the solution to \( Ax=b \) within the given tolerance tol. Otherwise, the last computed approximation to the solution is returned.
[in]	A	The square non-singular system matrix \( A \).
[in]	b	The right-hand side vector \( b \).

If the size of \( A \) is \( n \times n \), then the sizes of x and b must be \( n \) also. The vector x must have capacity \( n \).

Mandatory inputs to the BiCGstab algorithm:

Parameters

[in]	max_iterations	The maximum number of iterations this algorithm may perform.
[in]	tol	The relative tolerance which determines when an an approximated solution \( x \) becomes acceptable. Must be positive and non-zero.

Additional outputs of this algorithm:

Parameters

[out]	iterations	When grb::SUCCESS is returned, the number of iterations that were required to obtain an acceptable approximate solution.
[out]	residual	When grb::SUCCESS is returned, the square of the 2-norm of the residual; i.e., \( (r,r) \), where \( r = b - Ax \).

To operate, this algorithm requires a workspace consisting of six vectors of length and capacity \( n \). If vectors with less capacity are passed as arguments, grb::ILLEGAL will be returned.

Parameters

[in] r,rhat,p,v,s,t Workspace vectors required for BiCGstab.

The BiCGstab operates over a field defined by the following algebraic structures:

Parameters

[in]	semiring	Defines the domains as well as the additive and the multicative monoid.
[in]	minus	The inverse of the additive operator.
[in]	divide	The inverse of the multiplicative operator.

Note: When compiling with the _DEBUG macro defined, the print-out statements require sqrt as an additional algebraic concept. This concept presently lives "outside" of ALP.

Valid descriptors to this algorithm are:

descriptors::no_casting
descriptors::transpose

Returns: grb::SUCCESS If an acceptable solution is returned.; grb::FAILED If the algorithm failed to find an acceptable solution and returns an approximate one with the given residual.; grb::MISMATCH If two or more of the input arguments have incompatible sizes.; grb::MISMATCH If one or more of the workspace vectors has an incompatible size.; grb::ILLEGAL If tol is zero or negative.; grb::ILLEGAL If x has capacity less than \( n \).; grb::ILLEGAL If one or more of the workspace vectors has a capacity less than \( n \).; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ conjugate_gradient()

grb::RC grb::algorithms::conjugate_gradient	(	grb::Vector< IOType, backend > &	x,
		const grb::Matrix< NonzeroType, backend, RSI, RSI, NZI > &	A,
		const grb::Vector< InputType, backend > &	b,
		const size_t	max_iterations,
		ResidualType	tol,
		size_t &	iterations,
		ResidualType &	residual,
		grb::Vector< IOType, backend > &	r,
		grb::Vector< IOType, backend > &	u,
		grb::Vector< IOType, backend > &	temp,
		const Ring &	ring = `Ring()`,
		const Minus &	minus = `Minus()`,
		const Divide &	divide = `Divide()`
	)

Solves a linear system \( b = Ax \) with \( x \) unknown by the Conjugate Gradients (CG) method on general fields.

Does not perform any preconditioning– for preconditioned CG and full documentation, please see preconditioned_conjugate_gradient.

◆ cosine_similarity()

RC grb::algorithms::cosine_similarity	(	OutputType &	similarity,
		const Vector< InputType1 > &	x,
		const Vector< InputType2 > &	y,
		const Ring &	ring = `Ring()`,
		const Division &	div = `Division()`
	)

Computes the cosine similarity.

Given two vectors \( x, y \) of equal size \( n \), this function computes \( \alpha = \frac{ (x,y) }{ ||x||_2\cdot||y||_2 } \).

The 2-norms and inner products are computed according to the given semi- ring. However, the norms make use of the standard sqrt and so the algorithm assumes a regular field is used. Effectively, hence, the semiring controls the precision / data types under which the computation is performed.

Template Parameters

descr	The descriptor under which to perform the computation.
OutputType	The type of the output element (scalar).
InputType1	The type of the first vector.
InputType2	The type of the second vector.
Ring	The semiring used.
Division	Which binary operator correspond to division corresponding to the given Ring.

Parameters

[out]	similarity	Where to fold the result into.
[in]	x	The non-zero left-hand input vector.
[in]	y	The non-zero right-hand input vector.
[in]	ring	The semiring to compute over.
[in]	div	The division operator corresponding to ring.

Note: The vectors x and/or y may be sparse or dense.

The argument div is optional. It will map to grb::operators::divide by default.

Returns: grb::SUCCESS If the computation was successful.; grb::MISMATCH If the vector sizes do not match. The output similarity is untouched – the call to this algorithm will have no other effects than returning grb::MISMATCH.; grb::ILLEGAL In case x is all zero, and/or when y is all zero. The output similarity is undefined.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ gmres()

grb::RC grb::algorithms::gmres	(	grb::Vector< NonzeroType > &	x,
		const grb::Matrix< NonzeroType > &	A,
		const grb::Vector< NonzeroType > &	b,
		const size_t	n_restart,
		const size_t	max_iterations,
		const ResidualType	tol,
		size_t &	iterations,
		size_t &	iterations_gmres,
		size_t &	iterations_arnoldi,
		ResidualType &	residual,
		std::vector< grb::Vector< NonzeroType > > &	Q,
		std::vector< NonzeroType > &	Hmatrix,
		grb::Vector< NonzeroType > &	temp,
		std::vector< NonzeroType > &	temp3,
		const Ring &	ring = `Ring()`,
		const Minus &	minus = `Minus()`,
		const Divide &	divide = `Divide()`,
		const std::function< ResidualType(ResidualType) > &	sqrtX = `std_sqrt< ResidualType, ResidualType >`
	)

Solves a linear system \( b = Ax \) with \( x \) unknown using the Generalised Minimal Residual (GMRES) method on general fields.

Explicit preconditioning is possible by providing an initialised matrix \( M \) of matching size.

Template Parameters

descr	The user descriptor.
NonzeroType	The input/output vector/matrix nonzero type.
ResidualType	The type of the residual norm.
Ring	The semiring under which to perform GMRES.
Minus	The minus operator corresponding to the inverse of the additive operator of the given Ring.
Divide	The division operator corresponding to the inverse of the multiplicative operator of the given Ring.

Valid descriptors to this algorithm are:

descriptors::no_casting,
descriptors::transpose.

By default, i.e., if none of ring, minus, divide or sqrtX (nor their types) are explicitly provided by the user, the natural field on double data types will be assumed.

Note: An abstraction of a field that encapsulates Ring, Minus, and Divide may be more appropriate. This will also naturally ensure that demands on domain types are met.

GMRES algorithm inputs:

Parameters

[in,out]	x	On input: an initial guess to the solution.
[in]	A	The (square) indefinite system matrix.
[in]	b	The known right-hand side in \( Ax = b \). Must be structurally dense.

If A is \( n \times n \), then x and b must have matching length \( n \). The vector x furthermore must have a capacity of \( n \).

Additional GMRES algorithm inputs:

Parameters

[in]	max_iterations	The maximum number of GMRES iterations.
[in]	n_restart	The number of GMRES restart iterations, i.e., the maximum size of the Krylov subspace.
[in]	tol	The requested relative residual tolerance.

Additional outputs (besides x):

Parameters

[out]	iterations	Total number of interactions performed.
[out]	iterations_gmres	Number of GMRES interactions performed.
[out]	iterations_arnoldi	Number of Arnoldi interactions performed.
[out]	residual	Residual norm.

The GMRES algorithm requires four workspace buffers:

Parameters

[in,out]	Q	std::vector of length n_restart + 1. Each element of Q is grb::Vector of size \( n \). On output only fist n_restart vectors have meaning.
[in,out]	temp	A temporary vector of the same size as x.
[in,out]	Hmatrix	std::vector of length n_restart x n_restart used to store temporary data.
[in,out]	temp3	std::vector of length n_restart used to store temporary data.

Finally, the algebraic structures over which the GMRES is executed are given:

Parameters

[in]	ring	The semiring under which to perform the GMRES.
[in]	minus	The inverse of the additive operator of ring.
[in]	divide	The inverse of the multiplicative operator of ring.

Additional algebraic structure used by norm2 primitive:

Parameters

[in] sqrtX The square root (inverse of the square), not necessarily a closed operation on fields which describe vector norms. E.g., for complex field sqrtX maps real numbers to real numbers. If not provided explicitly the std::sqrt() is used is used if possible; if not, a compile time error is raised.

This algorithm may return one of the following error codes:

Returns: grb::SUCCESS When the algorithm has converged to a solution within the given max_iterations and tol.; grb::FAILED When the algorithm did not converge within the given max_iterations.; grb::ILLEGAL When A is not square.; grb::MISMATCH When x or b does not match the size of A.; grb::ILLEGAL When x does not have capacity \( n \).; grb::ILLEGAL When at least one of the workspace vectors does not have capacity \( n \).; grb::ILLEGAL If tol is not strictly positive.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ hessolve()

grb::RC grb::algorithms::hessolve	(	std::vector< NonzeroType > &	H,
		const DimensionType	n,
		const DimensionType &	kspspacesize,
		const ResidualType	tol,
		std::vector< NonzeroType > &	vecx,
		const Ring &	ring = `Ring()`,
		const Minus &	minus = `Minus()`,
		const Divide &	divide = `Divide()`,
		const std::function< ResidualType(ResidualType) > &	sqrtX = `std_sqrt< ResidualType, ResidualType >`
	)

Solves the least linear square problem of size n - 1, defined by the equation \( A x = b, \) where the \( A \) is an upper-Hessenberg matrix of size n x kspspacesize, and the vector \( b \), of length n, has only the first element nonzero.

The solution vector \( x \) is of length kspspacesize.

The algorithm first performs Givens rotations in order to reduce the upper Hessenberg matrix to upper rectangular form. After that, using back-substitution, the final solution \( x \) is calculated.

The matrix \( A \) and the vector \( b \) are stored in H, and the output vector \( x \) is copied into vecx.

Template Parameters

NonzeroType	The input/output vector/matrix nonzero type.
ResidualType	The type of the residual norm.
Ring	The semiring under which to perform GMRES.
Minus	The minus operator corresponding to the inverse of the additive operator of the given Ring.
Divide	The division operator corresponding to the inverse of the multiplicative operator of the given Ring.

Parameters

[in,out]	H	std::vector of length n x n used to store the vector \( b \) (first n elements of H ) and the matrix \( A \) ( remaining elements of H ) in a row-major orientation. On input only the first kspspacesize columns of \( A \) have meaning. On output, H is overwritten.
[out]	vecx	std::vector of length n that stores the solution \( x \). On input vecx is ignored. On output only the first kspspacesize elements have meaning. It is used to update the GMRES solution vector.

Note: The algorithm can do an in-place update of the vector \( b \) without the need to copy the solution into vecx.

Parameters

[in]	tol	The requested relative residual tolerance. Must be strictly positive.
[in]	n	The dimension of H and vecx, and the maximum size of the Krylov subspace.
[in]	kspspacesize	The size of Krylov subspace.

The algebraic structures over which the GMRES is executed:

Parameters

[in]	ring	The semiring under which to perform the GMRES.
[in]	minus	The inverse of the additive operator of ring.
[in]	divide	The inverse of the multiplicative operator of ring.

Additional algebraic structure used by norm2 primitive:

Parameters

[in] sqrtX The square root (inverse of the square), not necessarily a closed operation on fields which describe vector norms. E.g., for a complex field, sqrtX maps real numbers to real numbers. If not provided explicitly, the std::sqrt() is used if possible, and if not, a compile time error is raised.

This algorithm may return one of the following error codes:

Returns: grb::SUCCESS When the algorithm has solved the least linear square problem successfully.; grb::ILLEGAL When the size of H is smaller than \( n \times n \).; grb::ILLEGAL When kspspacesize is not in the interval \( \left[ 1, n \right> \).; grb::ILLEGAL If tol is not strictly positive.; grb::MISMATCH When vecx size in is smaller or equal to \( n \).; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS then becomes undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ kcore_decomposition()

RC grb::algorithms::kcore_decomposition	(	const Matrix< NZType > &	A,
		Vector< IOType > &	core,
		Vector< IOType > &	distances,
		Vector< IOType > &	temp,
		Vector< IOType > &	update,
		Vector< bool > &	status,
		IOType &	k
	)

The \( k \)-core decomposition algorithm.

Note: This algorithm is smoke-tested using a ground-truth output coreness vector corresponding to the EPA matrix. However, the ground truth was generated using an earlier version of this algorithm, run using an earlier version of ALP/GraphBLAS. This solution was manually verified against an external algorithm. A better testing methodology compares against a ground truth generated by such an external baseline– see GitHub issue #160, to which contributions would be warmly received.

Divides the input matrix into subgraphs with a coreness level. The coreness level \( k \) is defined as the largest subgraph in which each node has at least \( k \) neighbors in the subgraph.

Template Parameters

IOType	The value type of the \( k \)-core vectors, usually an integer type.
NZType	The type of the nonzero elements in the matrix.

Parameters

[in]	A	Matrix representing a graph with nonzero value at \( (i, j) \) an edge between node \( i \) and \( j \).
[out]	core	Empty vector of size and capacity \( n \). On output, if grb::SUCCESS is returned, stores the coreness level for each node.
[out]	k	The number of coreness lever that was found in the graph.

To operate, this algorithm requires a workspace of four vectors. The size and capacities of these must equal \( n \). The contents on input are ignored, and the contents on output are undefined. The work space consists of the buffer vectors distances, temp, update, and status.

Parameters

[in,out]	distances	Distance buffer
[in,out]	temp	First node update buffer
[in,out]	update	Second node update buffer
[in,out]	status	Finished/unfinished buffer

Returns: grb::SUCCESS If the coreness for all nodes are found.; grb::ILLEGAL If A is not square. All outputs are left untouched.; grb::MISMATCH If the dimensions of core or any of the buffer vectors does not match A. All outputs are left untouched.; grb::ILLEGAL If the capacity of one or more of core and the buffer vectors is less than \( n \).; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

If any non grb::SUCCESS error code is returned, then the contents of core are undefined, while k will be untouched by the algorithm.

Note: For undirected, unweighted graphs, use pattern matrix for A; i.e., use NZtype void; For unweighted graphs, IOType should be a form of unsigned integer. The value of any IOType element will be no more than the maximum degree found in the graph A.

Template Parameters

criticalSection The original MR had an eWiseLambda-based implementation that contains a critical section. This may or may not be faster than a pure ALP/GraphBLAS implementation, depending also on which backend is selected. Setting this template argument true selects the original eWiseLambda-based implementation, while otherwise a pure ALP/GraphBLAS implementation takes effect.

Note: In some non-exhaustive experiments, setting criticalSection to false leads to better performance on shared-memory parallel systems (using grb::reference_omp).

Warning: Setting criticalSection to true is not supported for the distributed-memory backends grb::BSP1D and grb::hybrid; see the corresponding code comment in the below algorithm for details.

For the above considerations, the default for criticalSection is presently set to false.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For additional performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

This algorithm is modelled after Li et al., "The K-Core Decomposition Algorithm Under the Framework of GraphBLAS", 2021 IEEE High Performance Extreme Computing Conference (HPEC), doi: 10.1109/HPEC49654.2021.9622845.

◆ kmeans_iteration()

RC grb::algorithms::kmeans_iteration	(	Matrix< IOType > &	K,
		Vector< std::pair< size_t, IOType > > &	clusters_and_distances,
		const Matrix< IOType > &	X,
		const size_t	max_iter = `1000`,
		const Operator &	dist_op = `Operator()`
	)

The kmeans iteration given an initialisation.

Parameters

[in,out]	K	k by m matrix containing the current k means as row vectors
[in]	clusters_and_distances	Vector containing the class and distance to centroid for each point
[in]	X	m by n matrix containing the n points to be classified as column vectors
[in]	max_iter	Maximum number of iterations
[in]	dist_op	Coordinatewise distance operator, squared difference by default

◆ knn()

RC grb::algorithms::knn	(	Vector< OutputType > &	u,
		const Matrix< InputType > &	A,
		const size_t	source,
		const size_t	k,
		Vector< bool > &	buf1
	)

Given a graph and a source vertex, indicates which vertices are contained within k hops.

This implementation is based on the matrix powers kernel over a Boolean semiring.

Parameters

[out]	u	The distance-k neighbourhood. Any prior contents will be ignored.
[in]	A	The input graph in (square) matrix form
[in]	source	The source vertex index.
[in]	k	The neighbourhood distance, or the maximum number of hops in a breadth-first search.

This algorithm requires the following workspace:

Parameters

[in,out] buf1 A buffer vector. Must match the size of A.

For \( n \times n \) matrices A, the capacity of u, buf1, and buf2 must equal \( n \).

Returns: grb::SUCCESS When the computation completes successfully.; grb::MISMATCH When the dimensions of u do not match that of A.; grb::MISMATCH If source is not in range of A.; grb::ILLEGAL If one or more of u, buf1, or buf2 has insufficient capacity.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ kpp_initialisation()

RC grb::algorithms::kpp_initialisation	(	Matrix< IOType > &	K,
		const Matrix< IOType > &	X,
		const Operator &	dist_op = `Operator()`
	)

a simple implementation of the k++ initialisation algorithm for kmeans

Parameters

[in,out]	K	k by m matrix containing the current k means as row vectors
[in]	X	m by n matrix containing the n points to be classified as column vectors
[in]	dist_op	Coordinatewise distance operator, squared difference by default

◆ label()

RC grb::algorithms::label	(	Vector< IOType > &	out,
		const Vector< IOType > &	y,
		const Matrix< IOType > &	W,
		const size_t	n,
		const size_t	l,
		const size_t	maxIterations = `1000`
	)

The label propagation algorithm.

Template Parameters

IOType The value type of the label vector. This will determine the precision of all computations this algorithm performs.

Parameters

[out]	out	The resulting labelled vector representing the n vertices.
[in]	y	Vector holding the initial labels from a total set of n vertices. The initial labels are assumed to correspond to the vertices corresponding to the first l entries of this vector. The labels must be either 0 ar 1.
[in]	W	Sparse symmetric matrix of size n by n, holding the weights between the n vertices. The weights must be positive (larger than 0). The matrix may be defective while the corresponding graph may not be connected.
[in]	n	The total number of vertices. If zero, and all of out, y, and W are empty, calling this function is equivalent to a no-op.
[in]	l	The number of vertices with an initial label. Must be larger than zero.
[in]	maxIterations	The maximum number of iterations this algorithm may execute. Optional. Default value: 1000.

Note: If the underlying graph is not connected then some components may be rendered immutable by this algorithm.

Returns: grb::ILLEGAL If one of the arguments passed to this function is illegal. The output will be left unmodified.; grb::ILLEGAL The vector out did not have full capacity available. The output will be left unmodified.; grb::ILLEGAL If n was nonzero but l was zero.; grb::SUCCESS If the computation converged within max iterations, or if n was zero.; grb::FAILED If the method did not converge within max iterations. The output will contain the latest iterand.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

[1] Kamvar, Haveliwala, Manning, Golub; ‘Extrapolation methods for accelerating the PageRank computation’, ACM Press, 2003.

◆ mpv()

RC grb::algorithms::mpv	(	Vector< IOType > &	u,
		const Matrix< InputType > &	A,
		const size_t	k,
		const Vector< IOType > &	v,
		Vector< IOType > &	temp,
		const Ring &	ring
	)

The matrix powers kernel.

Calculates \( y = A^k x \) for some integer \( k \geq 0 \) using the given semiring.

Template Parameters

descr	The descriptor used to perform this operation.
Ring	The semiring used.
IOType	The output vector type.
InputType	The nonzero type of matrix elements.
implementation	Which implementation to use.

Parameters

[out]	u	The output vector. Contents shall be overwritten. The supplied vector must match the row dimension size of A.
[in]	A	The square input matrix A. The supplied matrix must match the dimensions of u and v.
[in]	k	How many matrix–vector multiplications are requested.
[in]	v	The input vector v. The supplied vector must match the column dimension size of A. It may not be the same vector as u.
[in]	ring	The semiring to be used. This defines the additive and multiplicative monoids to be used.

This algorithm requires workspace:

Parameters

[in,out] temp A workspace buffer of matching size to the row dimension of A. May be the same vector as v, though note that the contents of temp on output are undefined.

This algorithm assumes that u and temp have full capacity. If this assumption does not hold, then a two-stage mpv must be employed instead.

Returns: grb::SUCCESS If the computation completed successfully.; grb::ILLEGAL If A is not square.; grb::MISMATCH If one or more of u, v, or temp has an incompatible size with A.; grb::ILLEGAL If one or more of u or temp does not have a full capacity.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.; grb::OVERLAP If one or more of v or temp is the same vector as u.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ norm2()

RC grb::algorithms::norm2	(	OutputType &	x,
		const Vector< InputType, backend, Coords > &	y,
		const Ring &	ring = `Ring()`,
		const std::function< OutputType(OutputType) >	sqrtX = `std_sqrt< OutputType, OutputType >`,
		const typename std::enable_if< std::is_floating_point< OutputType >::value, void >::type *	= `nullptr`
	)

Provides a generic implementation of the 2-norm computation.

Proceeds by computing a dot-product on itself and then taking the square root of the result.

This function is only available when the output type is floating point.

For return codes, exception behaviour, performance semantics, template and non-template arguments,

See also: grb::dot.

Parameters

[in,out]	x	On successful execution of this algorithm, on output, the 2-norm of y will have been added to x.
[in]	y	The vector to compute the norm of.
[in]	ring	The Semiring under which the 2-norm is to be computed.
[in]	sqrtX	The square root function which operates on real data type, as both input and output. If not explicitly provided, the std::sqrt() is used.

◆ preconditioned_conjugate_gradient()

grb::RC grb::algorithms::preconditioned_conjugate_gradient	(	grb::Vector< IOType, backend > &	x,
		const grb::Matrix< NonzeroType, backend, RSI, RSI, NZI > &	A,
		const grb::Vector< InputType, backend > &	b,
		const std::function< grb::RC(grb::Vector< IOType, backend > &, const grb::Vector< IOType, backend > &) > &	Minv,
		const size_t	max_iterations,
		ResidualType	tol,
		size_t &	iterations,
		ResidualType &	residual,
		grb::Vector< IOType, backend > &	r,
		grb::Vector< IOType, backend > &	u,
		grb::Vector< IOType, backend > &	temp,
		grb::Vector< IOType, backend > &	temp_precond,
		const Ring &	ring = `Ring()`,
		const Minus &	minus = `Minus()`,
		const Divide &	divide = `Divide()`
	)

Solves a preconditioned linear system \( b = M^{-1}Ax \) with \( x \) unknown by the Conjugate Gradients (CG) method on general fields.

The preconditioner M, alike to the system matrix A, must be symmetric positive definite.

Template Parameters

descr	The user descriptor
preconditioned	Whether to apply any given preconditioners.

Note: The default value for preconditioned is true and it is normally not necessary to override this value: if wishing to call a non-preconditioned CG, please see conjugate_gradient.

Template Parameters

IOType	The input/output vector nonzero type
ResidualType	The type of the residual
NonzeroType	The matrix nonzero type
InputType	The right-hand side vector nonzero type
Ring	The semiring under which to perform CG
Minus	The minus operator corresponding to the inverse of the additive operator of the given Ring.
Divide	The division operator corresponding to the inverse of the multiplicative operator of the given Ring.

Valid descriptors to this algorithm are:

descriptors::no_casting
descriptors::transpose

By default, i.e., if none of ring, minus, or divide (nor their types) are explicitly provided by the user, the natural field on double data types will be assumed.

Note: An abstraction of a field that encapsulates Ring, Minus, and Divide may be more appropriate. This will also naturally ensure that demands on domain types are met.

Parameters

[in,out]	x	On input: an initial guess to the solution. On output: the last computed approximation.
[in]	A	The (square) positive semi-definite system matrix.
[in]	b	The known right-hand side in \( Ax = b \). Must be structurally dense.

The preconditioner action \( M^{-1}r \) is supplied as a standard std::function that outputs a grb::RC and takes two arguments:

the output vector of the preconditioning action as a reference to a grb::Vector with element type IOType, and
the input vector to the preconditioner as a const-reference to a grb::Vector with the same element type. Both vectors are guaranteed to be dense when given to Minv.

A good preconditioner action Minv is both efficient to apply as well as drastrically reduces the number of CG iterations required. The latter is achieved by constructing Minv so that the condition number of \( M^{-1}A \) is much smaller than that of \( A \).

Parameters

[in] Minv The preconditioner action if preconditioned equals true.

If A is \( n \times n \), then x and b must have matching length \( n \). The vector x furthermore must have a capacity of \( n \).

CG algorithm inputs:

Parameters

[in]	max_iterations	The maximum number of CG iterations. Must be larger than zero.
[in]	tol	The requested relative tolerance.

Additional outputs (besides x):

Parameters

[out]	iterations	The number of iterations the algorithm has started.
[out]	residual	The residual corresponding to output x.

The CG algorithm requires three workspace buffers with capacity \( n \):

Parameters

[in,out]	r	A temporary vector of the same size as x.
[in,out]	u	A temporary vector of the same size as x.
[in,out]	temp	A temporary vector of the same size as x.
[in,out]	temp_precond	A temporary vector of the same size as x.

Note: If preconditioned is false, then both Minv and temp_precond are ignored. In this case, temp_precond need not have the same length as x, nor need it have full capacity.

Finally, the algebraic structures over which the CG is executed are given:

Parameters

[in]	ring	The semiring under which to perform the CG.
[in]	minus	The inverse of the additive operator of ring.
[in]	divide	The inverse of the multiplicative operator of ring.

This algorithm may return one of the following error codes:

Returns: grb::SUCCESS When the algorithm has converged to a solution within the given max_iterations and tol.; grb::FAILED When the algorithm did not converge within the given max_iterations.; grb::ILLEGAL When A is not square.; grb::MISMATCH When x or b does not match the size of A.; grb::ILLEGAL When x does not have capacity \( n \).; grb::ILLEGAL When at least one of the required workspace vectors does not have capacity \( n \).; grb::ILLEGAL If tol is not strictly positive.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

On output, the contents of the workspace r, u, and temp are always undefined. For non-grb::SUCCESS error codes, additional containers or states may be left undefined:

when grb::PANIC is returned, the entire program state, including the contents of all containers, become undefined;
when grb::ILLEGAL or grb::MISMATCH are returned and iterations equals zero, then all outputs are left unmodified compared to their contents at function entry;
when grb::ILLEGAL or grb::MISMATCH are returned and iterations is nonzero, then the contents of x are undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ preconditioned_gmres()

grb::RC grb::algorithms::preconditioned_gmres	(	grb::Vector< NonzeroType > &	x,
		const grb::Matrix< NonzeroType > &	M,
		const grb::Matrix< NonzeroType > &	A,
		const grb::Vector< NonzeroType > &	b,
		const size_t	n_restart,
		const size_t	max_iterations,
		const ResidualType	tol,
		size_t &	iterations,
		size_t &	iterations_gmres,
		size_t &	iterations_arnoldi,
		ResidualType &	residual,
		std::vector< grb::Vector< NonzeroType > > &	Q,
		std::vector< NonzeroType > &	Hmatrix,
		grb::Vector< NonzeroType > &	temp,
		std::vector< NonzeroType > &	temp3,
		const Ring &	ring = `Ring()`,
		const Minus &	minus = `Minus()`,
		const Divide &	divide = `Divide()`,
		const std::function< ResidualType(ResidualType) > &	sqrtX = `std_sqrt< ResidualType, ResidualType >`
	)

Solves a linear system \( b = Ax \) with \( x \) unknown using the Generalised Minimal Residual (GMRES) method on general fields.

Explicit preconditioning is enabled, using an initialised matrix \( M \) of matching size.

Template Parameters

descr	The user descriptor.
NonzeroType	The input/output vector/matrix nonzero type.
ResidualType	The type of the residual norm.
Ring	The semiring under which to perform GMRES.
Minus	The minus operator corresponding to the inverse of the additive operator of the given Ring.
Divide	The division operator corresponding to the inverse of the multiplicative operator of the given Ring.

Valid descriptors to this algorithm are:

descriptors::no_casting,
descriptors::transpose.

By default, i.e., if none of ring, minus, divide or sqrtX (nor their types) are explicitly provided by the user, the natural field on double data types will be assumed.

Note: An abstraction of a field that encapsulates Ring, Minus, and Divide may be more appropriate. This will also naturally ensure that demands on domain types are met.

GMRES algorithm inputs:

Parameters

[in,out]	x	On input: an initial guess to the solution.
[in]	M	The (square) preconditioning matrix.
[in]	A	The (square) indefinite system matrix.
[in]	b	The known right-hand side in \( Ax = b \). Must be structurally dense.

If A is \( n \times n \), then x and b must have matching length \( n \). The vector x furthermore must have a capacity of \( n \). The matrix M must be an \( n \times n \) square matrix.

Additional GMRES algorithm inputs:

Parameters

[in]	max_iterations	The maximum number of GMRES iterations.
[in]	n_restart	The number of GMRES restart iterations, i.e., the maximum size of the Krylov subspace.
[in]	tol	The requested relative residual tolerance.

Additional outputs (besides x):

Parameters

[out]	iterations	Total number of interactions performed.
[out]	iterations_gmres	Number of GMRES interactions performed.
[out]	iterations_arnoldi	Number of Arnoldi interactions performed.
[out]	residual	The residual norm.

The GMRES algorithm requires four workspace buffers:

Parameters

[in,out]	Q	std::vector of length n_restart + 1. Each element of Q is grb::Vector of size \( n \). On output only fist n_restart vectors have meaning.
[in,out]	temp	A temporary vector of the same size as x.
[in,out]	Hmatrix	std::vector of length n_restart x n_restart used to store temporary data.
[in,out]	temp3	std::vector of length n_restart used to store temporary data.

Finally, the algebraic structures over which the GMRES is executed are given:

Parameters

[in]	ring	The semiring under which to perform the GMRES.
[in]	minus	The inverse of the additive operator of ring.
[in]	divide	The inverse of the multiplicative operator of ring.

Additional algebraic structure used by norm2 primitive:

Parameters

[in] sqrtX The square root (inverse of the square), not necessarily a closed operation on fields which describe vector norms. E.g., for complex field sqrtX maps real numbers to real numbers. If not provided explicitly the std::sqrt() is used if possible, if not, a compile time error is raised.

This algorithm may return one of the following error codes:

Returns: grb::SUCCESS When the algorithm has converged to a solution within the given max_iterations and tol.; grb::FAILED When the algorithm did not converge within the given max_iterations.; grb::ILLEGAL When A is not square.; grb::MISMATCH When x or b or M does not match the size of A.; grb::ILLEGAL When x does not have capacity \( n \).; grb::ILLEGAL When at least one of the workspace vectors does not have capacity \( n \).; grb::ILLEGAL If tol is not strictly positive.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ simple_pagerank()

RC grb::algorithms::simple_pagerank	(	Vector< IOType > &	pr,
		const Matrix< NonzeroT > &	L,
		Vector< IOType > &	pr_next,
		Vector< IOType > &	pr_nextnext,
		Vector< IOType > &	row_sum,
		const IOType	alpha = `0.85`,
		const IOType	conv = `0.0000001`,
		const size_t	max = `1000`,
		size_t *const	iterations = `nullptr`,
		double *const	quality = `nullptr`
	)

The canonical PageRank algorithm.

Template Parameters

descr	The descriptor under which to perform the computation.
IOType	The value type of the pagerank vector. This will determine the precision of all computations this algorithm performs.
NonzeroT	The type of the elements of the nonzero matrix.

Parameters

[in,out]	pr	Vector of size and capacity \( n \), where \( n \) is the vertex size of the input graph L. On input, the contents of this vector will be taken as the initial guess to the final result, but only if the vector is dense; if it is, all entries of the initial guess must be nonzero, while it is not, this algorithm will make an initial guess. On output, if grb::SUCCESS is returned, the PageRank vector corresponding to L.
[in]	L	The input graph as a square link matrix of size \( n \).

To operate, this algorithm requires a workspace of three vectors. The size and capacities of these must equal \( n \). The contents on input are ignored, and the contents on output are undefined.

This algorithm does not explicitly materialise the Google matrix \( G = \alpha L + (1-\alpha)ee^T \) over which the power iterations are exectuted.

Parameters

[in,out]	pr_next	Buffer for the PageRank algorithm.
[in,out]	pr_nextnext	Buffer for the PageRank algorithm.
[in,out]	row_sum	Buffer for the PageRank algorithm.

The PageRank algorithm holds the following optional parameters:

Parameters

[in]	alpha	The scaling factor. The default value is 0.85. This value must be smaller than 1, and larger than 0.
[in]	conv	If the difference between two successive iterations, in terms of its one-norm, is less than this value, then the PageRank vector is considered converged and this algorithm exits successfully. The default value is \( 10^{-8} \). If this value is set to zero, then the algorithm will continue until max iterations are reached. May not be negative.
[in]	max	The maximum number of power iterations. The default value is 1000. This value must be larger than 0.

The PageRank algorithm reports the following optional output:

Parameters

[out]	iterations	If not `nullptr`, the number of iterations the call to this algorithm took will be written to the location pointed to.
[out]	quality	If not `nullptr`, the last computed residual will be written to the location pointed to.

Returns: grb::SUCCESS If the computation converged within max iterations.; grb::ILLEGAL If L is not square. All outputs are left untouched.; grb::MISMATCH If the dimensions of pr and L do not match. All outputs are left untouched.; grb::ILLEGAL If an invalid value for alpha, conv, or max was given. All outputs are left untouched.; grb::ILLEGAL If the capacity of one or more of pr, pr_next, pr_nextnext, or row_sum is less than \( n \). All outputs are left untouched.; grb::FAILED If the PageRank method did not converge to within the given tolerance conv within the given maximum number of iterations max. The output pr is the last computed approximation, while iterations and quality are likewise set to max and the last computed residual, respectively.; grb::PANIC If an unrecoverable error has been encountered. The output as well as the state of ALP/GraphBLAS is undefined.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ sparse_nn_single_inference() [1/2]

grb::RC grb::algorithms::sparse_nn_single_inference	(	grb::Vector< IOType > &	out,
		const grb::Vector< IOType > &	in,
		const std::vector< grb::Matrix< WeightType > > &	layers,
		const std::vector< BiasType > &	biases,
		grb::Vector< IOType > &	temp,
		const ReluMonoid &	relu = `ReluMonoid()`,
		const Ring &	ring = `Ring()`
	)

Performs an inference step of a single data element through a Sparse Neural Network defined by num_layers sparse weight matrices and num_layers biases.

The initial single data element may be sparse also, such as common in Graph Neural Networks (GNNs).

Inference here is a repeated sequence of application of a sparse linear layer, addition of a bias factor, and the application of a ReLU.

We employ a linear algebraic formulation where the ReLU and the bias application are jointly applied via a max-operator.

This formalism follows closely the linear algebraic approach to the related IEEE/MIT GraphChallenge problem, such as, for example, described in

Combinatorial Tiling for Sparse Neural Networks F. Pawlowski, R. H. Bisseling, B. Uçar and A. N. Yzelman 2020 IEEE High Performance Extreme Computing (HPEC) Conference

Parameters

[out]	out	The result of inference through the neural network.
[in]	in	The input vector, may be sparse or dense.
[in]	layers	A collection of linear layers. Each layer is assumed to be square and of the equal size to one another.

This implies that all layers are \( n \times n \). The vectors in and out hence must be of length \( n \).

Commonly, as an input propagates through a network, the features become increasingly dense. Hence out is assumed to have full capacity in order to potentially store a fully dense activation vector.

Inference proceeds under a set of biases, one for each layer. Activation vectors are added a constant bias value prior to applying the given relu function. This function does not perform tresholding.

Parameters

[in] biases An array of num_layers bias factors.

Inference is done using a single buffer that is alternated with out:

Parameters

[in,out] temp A buffer of size and capacity \( n \).

Finally, optional arguments define the algebraic structures under which inference proceeds:

Parameters

[in]	relu	The non-linear ReLU function to apply element-wise.
[in]	ring	The semiring under which to perform the inference.

The default algebraic structures are standard relu (i.e., max), min for tresholding, and the real (semi-) ring.

Valid descriptors for this algorithm are:

descriptor::no_casting

Note: This algorithm applies the propagation through layers in-place. To facilitate this, only square layers are allowed. Non-square layers would require the use of different vectors at every layer.

Returns: grb::SUCCESS If the inference was successful; grb::ILLEGAL If the size of layers does not match that of baises.; grb::MISMATCH If at least one pair of dimensions between layers, in, out, and temp do not match.; grb::ILLEGAL If at least one layer was not square.; grb::ILLEGAL If the capacities of one or more of out and temp were not full.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ sparse_nn_single_inference() [2/2]

grb::RC grb::algorithms::sparse_nn_single_inference	(	grb::Vector< IOType > &	out,
		const grb::Vector< IOType > &	in,
		const std::vector< grb::Matrix< WeightType > > &	layers,
		const std::vector< BiasType > &	biases,
		const ThresholdType	threshold,
		grb::Vector< IOType > &	temp,
		const ReluMonoid &	relu = `ReluMonoid()`,
		const MinMonoid &	min = `MinMonoid()`,
		const Ring &	ring = `Ring()`
	)

Performs an inference step of a single data element through a Sparse Neural Network defined by num_layers sparse weight matrices and num_layers biases.

The initial single data element may be sparse also, such as common in Graph Neural Networks (GNNs).

Inference here is a repeated sequence of application of a sparse linear layer, addition of a bias factor, and the application of a ReLU.

We employ a linear algebraic formulation where the ReLU and the bias application are jointly applied via a max-operator.

This formalism follows closely the linear algebraic approach to the related IEEE/MIT GraphChallenge problem, such as for example described in

Combinatorial Tiling for Sparse Neural Networks F. Pawlowski, R. H. Bisseling, B. Uçar and A. N. Yzelman 2020 IEEE High Performance Extreme Computing (HPEC) Conference

Parameters

[out]	out	The result of inference through the neural network.
[in]	in	The input vector, may be sparse or dense.
[in]	layers	A collection of linear layers. Each layer is assumed to be square and of the equal size to one another.

This implies that all layers are \( n \times n \). The vectors in and out hence must be of length \( n \).

Commonly, as an input propagates through a network, the features become increasingly dense. Hence out is assumed to have full capacity in order to potentially store a fully dense activation vector.

Inference proceeds under a set of biases, one for each layer. Activation vectors are added a constant bias value prior to applying the given relu function. After application, the resulting vector is furthermore tresholded. The treshold is assumed constant over all layers.

Parameters

[in]	biases	An array of num_layers bias factors.
[in]	threshold	The value used for thresholding.

Inference is done using a single buffer that is alternated with out:

Parameters

[in,out] temp A buffer of size and capacity \( n \).

Finally, optional arguments define the algebraic structures under which inference proceeds:

Parameters

[in]	relu	The non-linear ReLU function to apply element-wise.
[in]	min	Operator used for thresholding. Maximum feature value is hard-coded to 32, as per the GraphChallenge.
[in]	ring	The semiring under which to perform the inference.

The default algebraic structures are standard relu (i.e., max), min for tresholding, and the real (semi-) ring.

Valid descriptors for this algorithm are:

descriptor::no_casting

Note: This algorithm applies the propagation through layers in-place. To facilitate this, only square layers are allowed. Non-square layers would require the use of different vectors at every layer.; Thresholding here means that feature maps as propagated through the neural network are capped at some maximum value, threshold.

Returns: grb::SUCCESS If the inference was successful; grb::ILLEGAL If the size of layers does not match that of baises.; grb::MISMATCH If at least one pair of dimensions between layers, in, out, and temp do not match.; grb::ILLEGAL If at least one layer was not square.; grb::ILLEGAL If the capacities of one or more of out and temp were not full.

Performance semantics

This function does not allocate nor free dynamic memory, nor shall it make any system calls.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ spy() [1/3]

RC grb::algorithms::spy	(	grb::Matrix< IOType > &	out,
		const grb::Matrix< InputType > &	in
	)

Given an input matrix and a smaller output matrix, map nonzeroes from the input matrix into the smaller one and count the number of nonzeroes that are mapped from the bigger matrix into the smaller.

Template Parameters

normalize If set to true, will not compute a number of mapped nonzeroes, but its inverse instead (one divided by the count). The default value for this template parameter is false.

Parameters

[out]	out	The smaller output matrix.
[in]	in	The larger input matrix.

Returns: SUCCESS If the computation completes successfully.; ILLEGAL If out has a number of rows or columns larger than that of in.

Warning: Explicit zeroes (that when cast from InputType to bool read false) will be counted as a nonzero by this algorithm.

Note: To not count explicit zeroes, pre-process the input matrix in, for example, as follows: grb::set( tmp, in, true );, with tmp a Boolean or pattern matrix of the same size as in.

Performance semantics

Warning: This algorithm does NOT request workspace buffers since due to the use of level-3 primitives it will have to allocate anyway– as such, this algorithm does not have clear performance semantics and should be used with care.

For performance semantics regarding work, inter-process data movement, intra-process data movement, synchronisations, and memory use, please see the specification of the ALP primitives this function relies on. These performance semantics, with the exception of getters such as grb::nnz, are specific to the backend selected during compilation.

◆ spy() [2/3]

RC grb::algorithms::spy	(	grb::Matrix< IOType > &	out,
		const grb::Matrix< bool > &	in
	)

Specialisation for boolean input matrices in.

See grb::algorithms::spy.

◆ spy() [3/3]

RC grb::algorithms::spy	(	grb::Matrix< IOType > &	out,
		const grb::Matrix< void > &	in
	)

Specialisation for void input matrices in.

See grb::algorithms::spy.

◆ std_sqrt()

OutputType grb::algorithms::std_sqrt ( const InputType x )

An alias of std::sqrt where the input and output types are templated separately.

Template Parameters

OutputType	The output type of the square-root operation.
InputType	The input type of the square-root operation.

Parameters

[in] x The value to take the square root of.

Returns: The square root of x, cast to OutputType.

Relies on the standard std::sqrt implementation. If this is not available for InputType, the use of this operation will result in a compile-time error.

This operation is used as a default to the norm2 algorithm, as well as a default to algorithms that depend on it.

Namespaces

Classes

Functions

Detailed Description

Function Documentation

◆ bicgstab()

◆ conjugate_gradient()

◆ cosine_similarity()

◆ gmres()

◆ hessolve()

◆ kcore_decomposition()

◆ kmeans_iteration()

◆ knn()

◆ kpp_initialisation()

◆ label()

◆ mpv()

◆ norm2()

◆ preconditioned_conjugate_gradient()

◆ preconditioned_gmres()

◆ simple_pagerank()

◆ sparse_nn_single_inference() [1/2]

◆ sparse_nn_single_inference() [2/2]

◆ spy() [1/3]

◆ spy() [2/3]

◆ spy() [3/3]

◆ std_sqrt()