Julia drivers for PARDISO libraries. This interface is based on the type ParDiSO. This type contains the information needed by PARDISO to perform the essential operations.

Clone this repository in the julia packages folder:

cd $JULIA_PKGDIR                                        (default: ~/.julia/v0.4)
git clone https://github.com/verbof/PARDISO.jl PARDISO
cd PARDISO/src

Set correctly the FC variable to the correct value for your system in Makefile; then, modify the variable LIBPARDISO to poit at the correct PARDISO shared library needed for compiling.

The basic constructor

    ParDiSo(matrixtype::Integer, msglevel::Integer);

creates a variable that is going to be used for a matrix of the the type matrixtype and with a verbose-level msglevel (0: no messages from PARDISO, 1: verbose mode).

Since PARDISO is based on a CSR (Compressed Sparse Row) representation for sparse matrices, the Julia's SparseMatrixCSC{T<:Integer} <: AbstractSparseMatrix{T,Int32} have to be converted into SparsePardisoCSR format. In order to do so, the constructors

    SparsePardisoCSR{Float64}(A::SparseMatrixCSC{Float64,Int})
    SparsePardisoCSR{Complex128}(A::SparseMatrixCSC{Complex128,Int})

are provided. NOTE that this class saves just the upper triangular part of the original matrix if it belongs to one of the following classes:

• Real Symmetric Positive Definite
• Real Symmetric Positive Undefinite
• Complex Hermitian Positive Definite
• Complex Hermitian Undefinite
• Complex Symmetric

Before performing any operation using PARDISO, it is needed to initialize the ParDiSo object using the initPARDISO function.

Important: never change the values of pt.

The PARDISO solver relies on the matrix factorization to solve the systems. So, PARDISO needs to go through a preliminar factorization fase.

Checks if the input matrix is consistent according to the matrix-type.

    checkPARDISO{Tnzval}(pardiso::ParDiSo, A::SparsePardisoCSR{Tnzval})

Performs a reordering of the rows and a symbolic factorization on the matrix in order to reduce the fill-in in the factorization phase.

    smbfctPARDISO{Tnzval}(pardiso::ParDiSo, A::SparsePardisoCSR{Tnzval})

Computes an LU factorization for the matrix and stores internally the factors.

    factorPARDISO{Tnzval}(pardiso::ParDiSo, A::SparsePardisoCSR{Tnzval})

Solves the system Ax = b using the internal factors and performs an iterative refinement.

    solvePARDISO{Tnzval}(pardiso::ParDiSo, A::SparsePardisoCSR{Tnzval}, n_rhs::Int64, b::Vector{Tnzval})

Realeases internal PARDISO memory.

    freePARDISO(pardiso::ParDiSo)

using PARDISO;

M = 1000;
f = rand(M-1) + im*rand(M-1);
A = speye(M) + spdiagm(f,+1, M,M) + spdiagm(f,-1, M,M);      # Creates a COMPLEX SYMMETRIC matrix

pardiso = ParDiSo(+6, 1);

initPARDISO(pardiso);

X = rand(Complex128, M, 5);         # Random solution (5 Right-Hand Sides)
b = A*X;                            # RHS's
X = colwise(X);                     # column-wise representation (see PARDISO.jl)

A = SparsePardisoCSR(A);

println("PARDISO CHECK MATRIX");
checkPARDISO(pardiso, A);


println("SYMBOLIC FACTORIZATION + FACTORIZATION");
smbfctPARDISO(pardiso, A);
factorPARDISO(pardiso, A);

println();
println("SOLVE COMPLEX SYSTEM");
@time x = solvePARDISO(pardiso, A, b);                          # Solution computed by PARDISO

println("Residual:     ", norm(X-x)/norm(X));                   # Relative error on the solution
println("Total memory: ", memoryPARDISO(pardiso), " kylobites.\n");

freePARDISO(pardiso);