Direct Sparse Solver for Time Dependent Problems

[mandrew]

Hello,

I have a time-dependent problem that involves solving N quasi 9-diagonal matrices (top two rows full, with 9 nonzero diagonals) at each time step. I am currently using the direct sparse solver and I was wondering what the best method for calling DSS might be. I am storing all the nonzero elements for all N matrices in the same vector, then looping through the vector:

DO time=1,endtime
DO j=1,L

error = dss_factor_real(handle,MKL_DSS_DEFAULTS,values(((j-1)*(11*N-20)+1):(j*(11*N-20))) )
error = dss_solve_real(handle,MKL_DSS_DEFAULTS,rhs,nRhs,A(:,j))

END DO
END DO

Obviously, I can factor all the matrices prior to time-stepping, but I am not sure how this is stored in "handle." Can anyone give me some input on this issue? Also, if you know of a more efficient solver that would be great too.

Mandrew

[Sergey_P_Intel2]

Dear Mandrew,

First of all, DSS is the simplified interface to PARDISO, and it doesn't support many of PARDISO functionalities. Namely, PARDISO can work simultaneously with many matrices that have identical sparsity structure (see description of maxfct and mnum PARDISO parameters in MKL manual). Doing this way PARDISO keeps information about symbolic factorization and factors in one and the same handle. Consequently you can factorize all the matrices prior solving process. Note that only one reordering stage is required for the matrices.
You can also factorize all the matrices before their solving by using many handles (one per matrix). It is applicable to not only PARDISO but also DSS interface. In this way information about sparsity structure is duplicated in the handles, and you have to solve each matrix independently, applying reordering, factorization and solution phases.

With best regards,
Sergey

[mandrew]

Sergey,

Many thanks for the answer - I will read more about PARDISO.

Mandrew

[mandrew]

Quoting - Sergey Pudov (Intel)

You can also factorize all the matrices before their solving by using many handles (one per matrix). It is applicable to not only PARDISO but also DSS interface. In this way information about sparsity structure is duplicated in the handles, and you have to solve each matrix independently, applying reordering, factorization and solution phases.

How do I go about storing all the numerical factorizations of the matrices prior to beginning the simulation. In other words I would like to complete phase 1 and phase 2, since this will speed up the simulation. All I can figure is to do the following:

maxfct=L
phase=12

DO i=1,L

mnum=j

CALL pardiso(pt,maxfct,mnum,mtype,phase,N+1,values(((j-1)*(11*N-20)+1):(j*(11*N-20))),&
rowIndex,columns,idum,nRhs,iparm,msglvl,ddum,ddum,error)

END DO

Then when I begin the simulation I can do

phase=33

DO i=1,L

mnum=j

CALL pardiso(pt,maxfct,mnum,mtype,phase,N+1,values(((j-1)*(11*N-20)+1):(j*(11*N-20))),&
rowIndex,columns,idum,nRhs,iparm,msglvl,b,x,error)

END DO


Will this work, or is there something that I am not understanding?

Mandrew

[Sergey_P_Intel2]

Mandrew,

You are not quite right. As I described above, this PARDISO mode with many matrices that have one and the same sparse structure required only one reordering stage (phase=11): for the first matrix of the set. In your example, reordering stage (first part of phase=12) is called for all the matrices. So, please apply phase=12 to the first matrix in the set and phase=22 for the others.

Sergey

[Gennady_F_Intel]

Mandrew,

The problem you reported has been fixed into the latest MKL 10.2 Update 4 version released recently.

You can Download this version from Intel Registration Center: <>. See announcement about that here. Please let us know if the problem is still exist.

--Gennady