Civil engineer trying to program using MKL library Pardiso
I'm in the final stages of a finite element program, and need to solve the basic Ax=B matrix for x. I'm attempting to use MKL Pardiso, as the matrix is large and sparse. I've got the code operating up until the actual call to pardiso. Error I'm getting is "Unhandled exception at 0x000000014018b14b in doudican_current.exe: 0xC0000005: Access violation writing location 0x0000000050000163."Any thoughts? See code below. !*********************************************************************** ! SUBROUTINE DOUDICAN_SOLVE (GLK, BL, SOLUTION) ! ! USE all USE MKL_PARDISO ! implicit none ! REAL(KIND=DBL), INTENT(IN), DIMENSION(NEQ, NEQ) :: GLK REAL(KIND=DBL), INTENT(IN), DIMENSION(NEQ) :: BL REAL(KIND=DBL), INTENT(OUT), DIMENSION(NEQ) :: SOLUTION REAL(KIND=DBL), ALLOCATABLE, DIMENSION(:) :: Astar INTEGER, ALLOCATABLE, DIMENSION(:) :: JA INTEGER, DIMENSION(NEQ) :: IA INTEGER, DIMENSION(64) :: PT INTEGER, DIMENSION(64) :: IPARM INTEGER, DIMENSION(NEQ) :: PERM INTEGER :: MAXFCT, MNUM, MTYPE, PHASE, NRHS, MSGLVL, ERROR, M, N, & COUNT, COUNT1, COUNT2, COUNT2OLD, I, J ! ! MODIFY GLK TO BECOME Astar ! COUNT=0 DO I=1,NEQ DO J=1,NEQ IF (GLK(I,J)/=0) THEN COUNT=COUNT+1 END IF END DO END DO ! ALLOCATE ( Astar(COUNT), JA(COUNT) ) ! COUNT1=0 COUNT2=0 COUNT2OLD=0 DO I=1,NEQ COUNT2=COUNT2+1 DO J=1,NEQ IF (GLK(I,J)/=0) THEN COUNT1=COUNT1+1 Astar(COUNT1)=GLK(I,J) JA(COUNT1)=J IF (COUNT2OLD<COUNT2) THEN IA(COUNT2)=COUNT1 COUNT2OLD=COUNT2 END IF END IF END DO END DO ! ! DEFINE PARAMETERS ! DO I=1,64 PT(I)=0 IPARM(I)=0 END DO MAXFCT =1 MNUM = 1 MTYPE = 2 PHASE = 11 NRHS=1 MSGLVL = 0 ! ! CALL THE SPARSE MATRIX SOLVER PARDISO ! CALL PARDISO_ (PT, MAXFCT, MNUM, MTYPE, PHASE, NEQ, Astar, IA, JA, PERM, IPARM, MSGLVL, BL, SOLUTION, ERROR) ! RETURN END SUBROUTINE DOUDICAN_SOLVE ! !***********************************************************************
| |
Civil engineer trying to program using MKL library Pardiso
Recompile with /traceback and run again to correlate the IP with source file name and line number. We don't even know in which routine the crash occurred from the information given.
Show the declarations of the arguments to DOUDICAN_SOLVE that were made in the routine that called it.
| |
Civil engineer trying to program using MKL library Pardiso
You can see how little computer science experience I have. Please tell me if this is helpful or what you were looking for. LRC= 704 forrtl: severe (157): Program Exception - access violation Image PC Routine Line Source doudican_current. 000000014018B14B Unknown Unknown Unknown doudican_current. 000000014018AF31 Unknown Unknown Unknown doudican_current. 0000000140093B6F DOUDICAN_SOLVE 4599 doudican.f90 doudican_current. 0000000140007F95 MAIN__ 404 doudican.f90 doudican_current. 000000014051C09C Unknown Unknown Unknown doudican_current. 0000000140119BE7 Unknown Unknown Unknown doudican_current. 0000000140119AEE Unknown Unknown Unknown kernel32.dll 00000000773D652D Unknown Unknown Unknown ntdll.dll 0000000077ACC521 Unknown Unknown Unknown Line 4599 is the Call to Pardiso Line 404 is the Call from the main program to Doudican_Solve Declarations to the argumenst to Doudican_solve REAL(KIND=DBL), ALLOCATABLE, DIMENSION(:) :: BL, SOLUTION REAL(KIND=DBL), ALLOCATABLE, DIMENSION(:,:) :: GLK with KIND=DBL defined in a module as follows: MODULE all IMPLICIT NONE SAVE ! !Determine processor KIND number for single and double precision ! INTEGER, PARAMETER :: SGL = SELECTED_REAL_KIND(p=6,r=37) !p=decimal digits of precision INTEGER, PARAMETER :: DBL = SELECTED_REAL_KIND(p=13,r=200) !r=maximum range > 10^r
| |
Civil engineer trying to program using MKL library Pardiso
The traceback establishes that the crash was inside MKL/Pardiso. Usually, this is caused not by an error in Pardiso but by the user passing erroneous values of the index vectors IA, JA. Pardiso does check the input arguments, but it is possible to get erroneous input past those checks. The C0000005 trap is a second line of defence against such errors.
You have to examine whether the arguments to Pardiso satisfy the specifications, either using a debugger or by printing the contents just before the call.
I note at least one error: in the packed column storage convention, the start-of-column index array IA has dimension (NEQ+1). You have not filled in a value for the last element of the array. Here is a simpler version of the code for full matrix to packed column conversion:
KOUNT=COUNT(GLK /= 0d0)
!
ALLOCATE ( Astar(KOUNT), JA(KOUNT) )
!
KOUNT=0
IA(1)=1
DO J=1,NEQ ! process columns 1, 2, ...
IA(J+1)=IA(J) ! cumulate
DO I=1,NEQ
IF (GLK(I,J)/=0) THEN
IA(J+1)=IA(J+1)+1
KOUNT=KOUNT+1
Astar(KOUNT)=GLK(I,J)
JA(KOUNT)=I
END IF
END DO
END DO
How does subroutine Doudican_Solv know the value of NEQ?
| |
Civil engineer trying to program using MKL library Pardiso
also, first of all in such cases, we'd recommend to check the input matrix representation by setting parm(27)=1see the description of it into manual. --Gennady
| |
Civil engineer trying to program using MKL library Pardiso
Thanks mecej4 - good catch and fix on IA. I'm still debugging (error still hits), but to answer your question NEQ is passed in from Module "all". It is the dimension of the original sparse matrix (NEQ x NEQ).
| |
Civil engineer trying to program using MKL library Pardiso
Hi Gennady - thanks for the insight. Can I just enter the parameter definition as follows, or do I have to spec all 64 components of IPARM?! DEFINE PARAMETERS ! DO I=1,64 PT(I)=0 IPARM(I)=0 END DO IPARM(27)= 1 MAXFCT =1 MNUM = 1 MTYPE = 2 PHASE = 11 NRHS = 1 MSGLVL = 0
| |
Civil engineer trying to program using MKL library Pardiso
Thanks for all your help. After some more work, I've gotten it this far. The code below bonks out, giving a "Program Exception - Access Violation" at the Call to Pardiso. Any thoughts on why that might be? SUBROUTINE DOUDICAN_SOLVE (GLK, BL, SOLUTION) ! ! THIS SUBROUTINE ! USE all USE MKL_PARDISO ! implicit none ! REAL(KIND=DBL), INTENT(IN), DIMENSION(NEQ, NEQ) :: GLK REAL(KIND=DBL), INTENT(IN), DIMENSION(NEQ) :: BL REAL(KIND=DBL), INTENT(OUT), DIMENSION(NEQ) :: SOLUTION REAL(KIND=DBL), ALLOCATABLE, DIMENSION(:) :: Astar INTEGER, ALLOCATABLE, DIMENSION(:) :: JA INTEGER, DIMENSION(NEQ+1) :: IA !INTEGER, DIMENSION(64) :: PT TYPE(MKL_PARDISO_HANDLE), ALLOCATABLE :: PT(:) INTEGER, ALLOCATABLE :: IPARM(:) INTEGER, DIMENSION(NEQ) :: PERM INTEGER :: MAXFCT, MNUM, MTYPE, PHASE, NRHS, MSGLVL, ERROR, M, N, & COUNT, COUNT1, COUNT2, COUNT2OLD, I, J, KOUNT, DUMMY ! ! MODIFY GLK TO BECOME Astar ! KOUNT=COUNT(GLK/=0d0) ALLOCATE ( Astar(KOUNT), JA(KOUNT) ) KOUNT = 0 IA(1)=1 DO J=1,NEQ IA(J+1)=IA(J) DO I=1,NEQ IF (GLK(I,J)/=0) THEN IA(J+1)=IA(J+1)+1 KOUNT=KOUNT+1 Astar(KOUNT)=GLK(I,J) JA(KOUNT)=I END IF END DO END DO ! ! DEFINE PARAMETERS ! ALLOCATE ( PT (64) ) ALLOCATE ( IPARM (64) ) DO I=1,64 IPARM(I)=0 PT(I)%DUMMY=0 END DO ! IPARM(1) = 0 ! solver default IPARM(27) = 1 MTYPE = 2 MAXFCT =1 MNUM = 1 NRHS = 1 MSGLVL = 0 ERROR = 0 ! ! CALL THE SPARSE MATRIX SOLVER PARDISO ! PHASE = 13 ! CALL PARDISO_ (PT, MAXFCT, MNUM, MTYPE, PHASE, NEQ, Astar, IA, JA, PERM, NRHS, IPARM, MSGLVL, BL, SOLUTION, ERROR) ! ! RETURN END SUBROUTINE DOUDICAN_SOLVE
| | |