intel openmp和intel mkl冲突
环境:vs2015+ivf2017背景:程序采用mkl库pardiso并行进行矩阵计算,准备采用openmp对其他代码进行局部优化。
环境变量:NUMBER_OF_PROCESSORS=4,OMP_STACKSIZE=1000000
问题:单打开/Qmkl:parallel,程序正常运行
同时打开/Qopenmp_stubs和/Qmkl:parallel编译无问题,正常运行没有问题
同时打开/Qopenmp和/Qmkl:parallel编译无问题,在pardiso计算后,报stack损坏
现状:已搜索查询很多网站,还是没有发现如何规避此问题,求助大神解答。
报错
代码:
MODULE MMkl_pardiso
USE Constants
IMPLICIT NONE
SAVE
PRIVATE
PUBLIC:: cal_coomatrix_mkl_pardiso ! 非对称阵
CONTAINS
SUBROUTINE cal_coomatrix_mkl_pardiso(nnz, acoo, rowind, colind, m, n,b, x)
IMPLICIT NONE
INTEGER(KINT) :: job(8)
INTEGER(KINT) :: nnz
REAL(KREAL), INTENT(INOUT) :: acoo(:)
INTEGER(KINT),INTENT(INOUT) :: rowind(:)
INTEGER(KINT),INTENT(INOUT) :: colind(:)
INTEGER(KINT),INTENT(IN) :: m
INTEGER(KINT),INTENT(INOUT) :: n
REAL(KREAL), ALLOCATABLE :: acsr(:)
INTEGER(KINT),ALLOCATABLE :: ja(:)
INTEGER(KINT),ALLOCATABLE :: ia(:)
INTEGER(KINT) :: info
INTEGER(KINT) :: acsr_dim
!
REAL(KREAL), INTENT(INOUT) :: b(:)
REAL(KREAL), INTENT(INOUT) :: x(:)
!
!=====================the variables of pardiso====================
integer(kind=4) :: pt(64) ! kind=8:x64; kind=4:win32
integer :: maxfct, mnum, mtype, phase, nrhs, error, msglvl
integer, allocatable :: perm(:)
integer :: iparm(64)
integer, external :: mkl_get_max_threads
!=====================the variables of pardiso====================
job(1) = 1 ! if job(1)=1, the matrix in the coordinate format is converted to the CSR format.
job(2) = 1 ! if job(2)=1, one-based indexing for the matrix in CSR format is used.
job(3) = 1 ! if job(3)=1, one-based indexing for the matrix in coordinate format is used.
job(4) = 0
job(5) = nnz ! nzmax - maximum number of the non-zero elements allowed if job(1)=0.
job(6) = 0 ! all arrays acsr, ja, ia are filled in for the output storage.
acsr_dim = nnz
! n = nnz
ALLOCATE( acsr(acsr_dim) )
ALLOCATE( ja(acsr_dim) )
ALLOCATE( ia(m+1))
ALLOCATE( perm(n))
acsr = REAL_ZERO
ja = INT_ZERO
ia = INT_ZERO
CALL convert_coo_to_csr(job, n, acsr, ja, ia, nnz, acoo, rowind, colind, info)
!===========================solve equations=============================
mtype = 11! mtype = -2: real and symmetric indefinite,
! mtype =2: real and symmetric positive definite
! mtype = 11: real and nonsymmetric
! Initialize Intel MKL PARDISO with default parameters in accordance with the matrix type.
call pardisoinit (pt, mtype, iparm)
! pt: Solver internal data address pointer
maxfct = 1 ! Maximal number of factors in memory
! Generally used value is 1
mnum = 1 ! The number of matrix (from 1 to maxfct) to solve
! Generally used value is 1
! mtype =
! 1 Real and structurally symmetric
! 2 Real and symmetric positive definite
! -2 Real and symmetric indefinite
! 3 Complex and structurally symmetric
! 4 Complex and Hermitian positive definite
! -4 Complex and Hermitian indefinite
! 6 Complex and symmetric matrix
! 11 Real and nonsymmetric matrix
! 13 Complex and nonsymmetric matrix
! phase =
! 11 Analysis
! 12 Analysis, numerical factorization
! 13 Analysis, numerical factorization, solve
! 22 Numerical factorization
! 23 Numerical factorization, solve
! 33 Solve, iterative refinement
! 331 phase=33, but only forward substitution
! 332 phase=33, but only diagonal substitution
! 333 phase=33, but only backward substitution
! 0 Release internal memory for L and U of the matrix number mnum
! -1 Release all internal memory for all matrices
! n Number of equations in the sparse linear system A*X = B
! a(*) Contains the non-zero elements of the coefficient matrix A
! ia(n+1) rowIndex array in CSR3 format
! ja(*) columns array in CSR3 format
! perm(n) Holds the permutation vector of size n or specifies elements used for computing a partial solution
! nrhs Number of right-hand sides that need to be solved for
! iparm(64) This array is used to pass various parameters to Intel MKL PARDISO
! and to return some useful information after execution of the solver
! (see pardiso iparm Parameter for more details)
! msglvl 0 Intel MKL PARDISO generates no output
! 1 Intel MKL PARDISO prints statistical information
! b(n*nrhs) Right-hand side vectors
! x(n*nrhs) Solution vectors
! error Error indicator
! call pardiso (pt, maxfct, mnum, mtype, phase, n, a, ia, ja, perm, nrhs, iparm, msglvl, b, x, error)
perm = 0
nrhs = 1 ! Number of right-hand sides that need to be solved for.
x = 0.0D0
! iparm(1) = 0! iparm use default values
! iparm(3): parallel threads are obtained by MKL. In general, the number of cpus is the same as the number of threads
!call mkl_set_num_threads(4)
iparm(3) = mkl_get_max_threads()
error = 0
msglvl = 0
phase = 12! LU decompose, Analysis, numerical factorization
call pardiso( pt, maxfct, mnum, mtype, phase, m, acsr, ia, ja, perm, nrhs, iparm, msglvl, b, x, error )
phase = 33! Solve, iterative refinement
call pardiso( pt, maxfct, mnum, mtype, phase, m, acsr, ia, ja, perm, nrhs, iparm, msglvl, b, x, error )
phase = -1! Release all internal memory for all matrices
call pardiso (pt, maxfct, mnum, mtype, phase, n, acsr, ia, ja, perm, nrhs, iparm, msglvl, b, x, error)
DEALLOCATE( acsr )
DEALLOCATE( ja )
DEALLOCATE( ia)
DEALLOCATE( perm)
RETURN
END SUBROUTINE cal_coomatrix_mkl_pardiso
SUBROUTINEconvert_coo_to_csr(job, n, acsr, ja, ia, nnz, acoo, rowind, colind, info)
IMPLICIT NONE
INTEGER(KINT),INTENT(IN OUT):: job(8)
INTEGER(KINT),INTENT(IN OUT):: n
REAL(KREAL) :: acsr(:)
INTEGER(KINT) :: ja(:)
INTEGER(KINT) :: ia(:)
INTEGER(KINT) :: nnz
REAL(KREAL), INTENT(IN OUT):: acoo(:)
INTEGER(KINT),INTENT(IN OUT):: rowind(:)
INTEGER(KINT),INTENT(IN OUT):: colind(:)
INTEGER(KINT) :: info
call mkl_dcsrcoo(job, n, acsr, ja, ia, nnz, acoo, rowind, colind, info)
RETURN
END SUBROUTINE convert_coo_to_csr
END MODULE MMkl_pardiso
本帖最后由 li913 于 2019-11-1 15:30 编辑
你给的代码不全,没法执行。试试关闭嵌套并行。
integer(kind=4) :: pt(64) ! kind=8:x64; kind=4:win32
最后发现并非mkl和openmp冲突,是由于64位编译下pt需要用kind=8,问题解决
页:
[1]