Erratic behavior of cluster_sparse_solver

[OP1]

The code attached below will either crash (access violation error) or produce the right result, or an approximate result depending on the number of MPI processes it runs on. This code is run on a Win 7 64-bit machine with dual six-core processor.

The 8-by-8 matrix used in this example is Hermitian, positive definite.

One should probably note that the MKL documentation page for 'cluster_sparse_solver' (the definition of the argument "a", for instance) does not mention MTYPE=4; this omission is found in multiple other places as well (description of ipar(10) in the page 'cluster_sparse_solver iparm Parameter', for instance).

PROGRAM MAIN
USE MPI
USE MKL_CLUSTER_SPARSE_SOLVER
IMPLICIT NONE


INTEGER(KIND=4)  :: MAXFCT,MNUM,MTYPE,PHASE,N,IA(9),JA(13),PERM(1),NRHS,IPARM(64),MSGLVL,ERROR,RANK,I
COMPLEX(KIND=8)  :: A(13),B(8),X(8),RES,RES0,ALPHA,BETA,BS(8)
CHARACTER(LEN=1) :: MATDESCRA(6)
TYPE(MKL_CLUSTER_SPARSE_SOLVER_HANDLE) :: PT(64)

CALL MPI_INIT(ERROR)
CALL MPI_COMM_RANK(MPI_COMM_WORLD,RANK,ERROR)

MAXFCT = 1
MNUM   = 1
MTYPE  = 4

N      = 8
A(1)  = 15
A(2)  = (3,1)
A(3)  = 14
A(4)  = (2.0,1.1)
A(5)  = 2
A(6)  = (1.0,0.5)
A(7)  = 25
A(8)  = 20
A(9)  = 15
A(10) = (1.0,0.2)
A(11) = 17
A(12) = (2.0,0.5)
A(13) = 25
IA     = [1,3,5,7,8,9,11,13,14]
JA     = [1,2,2,3,3,4,4,5,6,7,7,8,8]

PERM   = 0
NRHS   = 1

MSGLVL = 1
DO I=1,64
    PT(I)%DUMMY = 0
ENDDO

IPARM     = 0
IPARM(1)  = 1   ! NO SOLVER DEFAULT
IPARM(2)  = 2   ! FILL-IN REORDERING FROM METIS
IPARM(6)  = 0   ! =0 SOLUTION ON THE FIRST N COMPOMENTS OF X
IPARM(8)  = 2   ! NUMBERS OF ITERATIVE REFINEMENT STEPS
IPARM(10) = 13  ! PERTURBE THE PIVOT ELEMENTS WITH 1E-13
IPARM(11) = 1   ! USE NONSYMMETRIC PERMUTATION AND SCALING MPS
IPARM(13) = 1   ! MAXIMUM WEIGHTED MATCHING ALGORITHM IS SWITCHED-OFF
IPARM(40) = 0   ! INPUT: MATRIX/RHS/SOLUTION STORED ON MASTER

PHASE  = 11
CALL CLUSTER_SPARSE_SOLVER(PT,MAXFCT,MNUM,MTYPE,PHASE,N,A,IA,JA,PERM,NRHS,IPARM,MSGLVL,B,X,MPI_COMM_WORLD,ERROR)
IF (RANK==0) WRITE(*,*) 'REORDERING COMPLETED ... '
IF (ERROR/=0) THEN
    IF (RANK==0) WRITE(*,*) 'THE FOLLOWING ERROR WAS DETECTED: ', ERROR
    CALL MPI_FINALIZE(ERROR)
    STOP 1
ENDIF

PHASE = 22
CALL CLUSTER_SPARSE_SOLVER(PT,MAXFCT,MNUM,MTYPE,PHASE,N,A,IA,JA,PERM,NRHS,IPARM,MSGLVL,B,X,MPI_COMM_WORLD,ERROR)
IF (RANK==0) WRITE(*,*) 'FACTORIZATION COMPLETED ... '
IF (ERROR/=0) THEN
    IF (RANK==0) WRITE(*,*) 'THE FOLLOWING ERROR WAS DETECTED: ',ERROR
    CALL MPI_FINALIZE(ERROR)
    STOP 1
ENDIF

PHASE = 33
IF (RANK==0) THEN
    DO I=1,N
        B(I) = (1,1)
    ENDDO
ENDIF
CALL CLUSTER_SPARSE_SOLVER(PT,MAXFCT,MNUM,MTYPE,PHASE,N,A,IA,JA,PERM,NRHS,IPARM,MSGLVL,B,X,MPI_COMM_WORLD,ERROR)
IF (ERROR/=0) THEN
    IF (RANK==0) WRITE(*,*) 'THE FOLLOWING ERROR WAS DETECTED: ',ERROR
    CALL MPI_FINALIZE(ERROR)
    STOP 1
ENDIF

MATDESCRA(1) ='H'
MATDESCRA(2) ='U'
MATDESCRA(3) ='N'
MATDESCRA(4) ='F'
ALPHA = 1.0
BETA  = 0.0
      
IF (RANK==0) THEN
    WRITE(*,*) 'SOLVE COMPLETED ... '
    WRITE(*,*) 'THE SOLUTION OF THE SYSTEM IS '
    DO I = 1, N
        WRITE(*,*) ' X(',I,') = ', X(I)
    ENDDO
    CALL MKL_ZCSRMV('N',N,N,ALPHA,MATDESCRA,A,JA,IA(1),IA(2),X,BETA,BS)
    RES  = (0.D0,0.D0)
    RES0 = (0.D0,0.D0)
    DO I=1,N
        RES = RES + (BS(I)-B(I))*CONJG((BS(I)-B(I)))
        RES0 = RES0 + B(I)*CONJG(B(I))
    ENDDO
    PRINT *, 'RELATIVE RESIDUAL = ', SQRT(ABS(RES))/SQRT(ABS(RES0))
ENDIF


END PROGRAM MAIN

The compile and link command is shown here:

Compiling with Intel(R) Visual Fortran Compiler 17.0 [Intel(R) 64]...
ifort /nologo /O2 /I"C:\Program Files (x86)\IntelSWTools\compilers_and_libraries_2017.0.048\windows\mkl\include" /I"C:\Program Files (x86)\IntelSWTools\mpi\2017.0.048\intel64\include" /module:"x64\Release\\" /object:"x64\Release\\" /Fd"x64\Release\vc120.pdb" /libs:dll /threads /c /Qlocation,link,"C:\Program Files (x86)\Microsoft Visual Studio 12.0\VC\\bin\amd64" /Qm64 "D:\TEMP\CLUSTER_SPARSE_SOLVER\MAIN.F90"
Linking...
Link /OUT:"x64\Release\CLUSTER_SPARSE_SOLVER.exe" /INCREMENTAL:NO /NOLOGO /LIBPATH:"C:\Program Files (x86)\IntelSWTools\compilers_and_libraries_2017.0.048\windows\mkl\lib\intel64_win" /LIBPATH:"C:\Program Files (x86)\IntelSWTools\compilers_and_libraries_2017.0.048\windows\mpi\intel64\lib" /MANIFEST /MANIFESTFILE:"x64\Release\CLUSTER_SPARSE_SOLVER.exe.intermediate.manifest" /MANIFESTUAC:"level='asInvoker' uiAccess='false'" /SUBSYSTEM:CONSOLE /IMPLIB:"D:\TEMP\CLUSTER_SPARSE_SOLVER\x64\Release\CLUSTER_SPARSE_SOLVER.lib" mkl_intel_lp64_dll.lib mkl_core_dll.lib mkl_intel_thread_dll.lib mkl_blacs_lp64_dll.lib impi.lib -qm64 "x64\Release\mkl_cluster_sparse_solver.obj" "x64\Release\MAIN.obj"
Embedding manifest...
mt.exe /nologo /outputresource:"D:\TEMP\Sandbox\CLUSTER_SPARSE_SOLVER\x64\Release\CLUSTER_SPARSE_SOLVER.exe;#1" /manifest "x64\Release\CLUSTER_SPARSE_SOLVER.exe.intermediate.manifest"

CLUSTER_SPARSE_SOLVER - 0 error(s), 0 warning(s)

(edited to change mkl_intel_lp64.lib to mkl_intel_lp64_dll.lib - this was a typo but does not change the outcome)

[Gennady_F_Intel]

Oliver, thanks for the report. We will check both of these problems. Have you had chance to check if the same problem will happen on Linux OS* too?

[OP1]

We are also having difficulties on a Linux system. For the program above the access violation occurs somewhere in MKL_CORE.DLL.

[OP1]

Also, the documentation for CLUSTER_SPARSE_SOLVER indicates that the first argument (PT) should be an array of 64 INTEGER(8) values. This does not match with the interface to the function, which requires PT to be a object of type MKL_CLUSTER_SPARSE_SOLVER_HANDLE.