SenK documentation

SenK is a C++ library for high-performance linear solvers that primarily provides solvers for sparse linear systems with CSR matrix storage format. You can easily apply SenK to your applicatioins because it is a header-only library that provides basic functions.

Requirement

I try to implement it in a way that it does not depend on the compiler and OS, but it has been tested on only macOS, Ubuntu, and CentOS. Since I am developing with C++17, please use a compiler that supports C++17. Installation of other libraries is not required because only C++ standard libraries are used.

Usage

Since SenK is a header-only library, no pre-build is required. Just include senk.hpp in your source file and specify the include path at compile time.
The SenK library can be downloaded from the GitHub Repository.

#include "senk.hpp"
int main(int argc, char **argv) {
    //Blhablhablha
}

g++ -O3 -std=c++17 -fopenmp main.cpp -I(PATH to SenK)/include -o a.out

Sample

A sample program to execute the ILU(0) preconditioned GMRES(m) solver in double-precision is as follows:

#include <string>
 
#include "senk.hpp"
#include "senk_test.hpp"
 
int main(int argc, char **argv)
{
    std::string filename = "PATH to a MatrixMarket file";
    // Read the MatrixMarket file.
    double *tval;
    int *tcind, *trptr;
    int N, M;
    senk::io::Shape shape;
    bool removeZeros = true;
    senk::io::ReadMatrixMarket(
        filename, &tval, &tcind, &trptr,
        &N, &M, &shape, removeZeros);
    
    // Expand the input matrix if it is symmetric.
    double *val;
    int *cind, *rptr;
    if(shape == senk::io::Shape::Sym) {
        senk::matrix::Expand(tval, tcind, trptr, &val, &cind, &rptr, N);
    }else {
        senk::matrix::Duplicate(tval, tcind, trptr, &val, &cind, &rptr, N);
    }
    if(!senk::matrix::CheckStructure(val, cind, rptr, N)) exit(1);
 
    double *b = senk::utils::SafeMalloc<double>(N);
    double *x = senk::utils::SafeMalloc<double>(N);
    senk::utils::Set<double>(1, b, N);
    senk::utils::Set<double>(0, x, N);
 
    // Perform the ILU(0) factorization.
    senk::matrix::Duplicate(val, cind, rptr, &tval, &tcind, &trptr, N);
    senk::matrix::Ilu0(tval, tcind, trptr, N);
 
    // Split the factorized matrix into L and U.
    double *lval;
    int *lcind, *lrptr;
    double *uval;
    int *ucind, *urptr;
    senk::matrix::Split(
        tval, tcind, trptr,
        &lval, &lcind, &lrptr,
        &uval, &ucind, &urptr,
        nullptr, N, "L-DU", true);
 
    // Execute the ILU(0)-GMRES(m) solver.
    int max_iter = 6000;
    int m = 50;
    double epsilon = 1.0e-8;
    double nrm_b = senk::blas1::Nrm2<double>(b, N);  
    senk::solver::IluGmresm<double>(
        val, cind, rptr,
        lval, lcind, lrptr, uval, ucind, urptr,
        b, x, nrm_b, max_iter/m, m, N, epsilon);
    
    senk::test::RelativeResidualError(
        val, cind, rptr, b, x, N);
    
    return 0;
}

Performance evaluations

Coming soon.

Implementation status

• Solver
  • Completed: GMRES(m), GCR(m), BiCGStab method
  • Ongoing: CG method
• Preconditioner
  • Completed: ILU(0), ILU(p), ILUB preconditioner
  • Ongoing: IC(0), IC(p) preconditioner

ILUB preconditioner

ILUB preconditioning is one of the preconditioning techniques based on ILU factorization. The control of fill-in during ILU factorization plays an important role. Since there is a trade-off between increased computational cost due to fill-ins and improved convergence property, it is necessary to control fill-in appropriately. In ILUB preconditioning, fill-in is controlled based on a certain (vertical) block structure. It thereby efficiently utilizes SIMD operations for calculations inside the block during forward/backward substitution in the preconditioning process. The SIMD operations can mitigate the effects of increased computational complexity, and the benefits of improved convergence property are more readily apparent.

Author

Kengo Suzuki

Date

5/9/2022