topic dss_solve_real in Intel® oneAPI Math Kernel Library

dss_solve_real

siva — Mon, 16 Mar 2009 21:51:12 GMT

Hi,

I am trying to do the following with the MKL direct sparse solver. (I am using Intel Fortran 11.0.066 and MKL 10.1.0)

I have a symmetricpositive definite sparse matrix. After factoring it, I do

1. solve with MKL_DSS_FORWARD_SOLVE optionfollowed by solve with MKL_DSS_BACKWARD_SOLVE option
2. solve with MKL_DSS_DEFAULTS options

They should both give the same answer, but they don't.

I have attached a file where the matrix I am using is diagonal A = [2 0; 0 2]

Could you please tell me what I am doing wrong?

Thank you.

Siva

[cpp]include 'mkl_dss.f90'

program main

use mkl_dss

implicit none

type(MKL_DSS_HANDLE) :: hDSS

double precision :: values(3)
integer :: rowIndex(3), columns(3)
integer :: perm(2)

double precision :: rhs(2), x1(2), x2(2), x(2)

values = 0.d0
values(1) = 2.d0
values(2) = 0.d0
values(3) = 2.d0

columns(1) = 1
columns(2) = 2
columns(3) = 2

rowIndex(1) = 1
rowIndex(2) = 3
rowIndex(3) = 4

rhs(1) = 1.d0
rhs(2) = 0.d0

perm(1) = 1
perm(2) = 2

if (dss_create( hDSS, MKL_DSS_DEFAULTS ) .ne. MKL_DSS_SUCCESS) then
    print*,'Error creating'
    stop
endif

if (dss_define_structure( hDSS, MKL_DSS_SYMMETRIC, rowIndex, 2, 2, columns, 3) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error defining structure'
    stop
endif    

if (dss_reorder( hDSS, MKL_DSS_AUTO_ORDER, perm ) .ne. MKL_DSS_SUCCESS)  then
    print*,'Error reordering'
    stop
endif

if (dss_factor_real( hDSS, MKL_DSS_POSITIVE_DEFINITE, values ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error factoring'
    stop
endif

if (dss_solve_real( hDSS, MKL_DSS_FORWARD_SOLVE, rhs, 1, x1 ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error forward solving'
    stop
endif

print*,x1

if (dss_solve_real( hDSS, MKL_DSS_BACKWARD_SOLVE, x1, 1, x2 ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error backward solving'
    stop
endif

print*,x2

if (dss_solve_real( hDSS, MKL_DSS_DEFAULTS, rhs, 1, x ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error solving'
    stop
endif

print*,x

end program main[/cpp]

Re: dss_solve_real

siva — Tue, 17 Mar 2009 02:30:10 GMT

Got it, have to do

MKL_DSS_FORWARD_SOLVE+MKL_DSS_REFINEMENT_OFF

etc. (see pardiso and dss documentation)

Quoting - siva@colorado.edu

[cpp]include 'mkl_dss.f90'

program main

use mkl_dss

implicit none

type(MKL_DSS_HANDLE) :: hDSS

double precision :: values(3)
integer :: rowIndex(3), columns(3)
integer :: perm(2)

double precision :: rhs(2), x1(2), x2(2), x(2)

values = 0.d0
values(1) = 2.d0
values(2) = 0.d0
values(3) = 2.d0

columns(1) = 1
columns(2) = 2
columns(3) = 2

rowIndex(1) = 1
rowIndex(2) = 3
rowIndex(3) = 4

rhs(1) = 1.d0
rhs(2) = 0.d0

perm(1) = 1
perm(2) = 2

if (dss_create( hDSS, MKL_DSS_DEFAULTS ) .ne. MKL_DSS_SUCCESS) then
    print*,'Error creating'
    stop
endif

if (dss_define_structure( hDSS, MKL_DSS_SYMMETRIC, rowIndex, 2, 2, columns, 3) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error defining structure'
    stop
endif    

if (dss_reorder( hDSS, MKL_DSS_AUTO_ORDER, perm ) .ne. MKL_DSS_SUCCESS)  then
    print*,'Error reordering'
    stop
endif

if (dss_factor_real( hDSS, MKL_DSS_POSITIVE_DEFINITE, values ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error factoring'
    stop
endif

if (dss_solve_real( hDSS, MKL_DSS_FORWARD_SOLVE, rhs, 1, x1 ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error forward solving'
    stop
endif

print*,x1

if (dss_solve_real( hDSS, MKL_DSS_BACKWARD_SOLVE, x1, 1, x2 ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error backward solving'
    stop
endif

print*,x2

if (dss_solve_real( hDSS, MKL_DSS_DEFAULTS, rhs, 1, x ) &
                                         .ne. MKL_DSS_SUCCESS)  then
    print*,'Error solving'
    stop
endif

print*,x

end program main[/cpp]

Re: dss_solve_real

Sergey_P_Intel2 — Tue, 17 Mar 2009 09:26:20 GMT

Hi, Siva!

Thank you for the question. Youhave found in MKL documentation the correct solution of the problem. Currently by default DSSperforms up to tworefinementstepsfor better accuracy.Butafter forward (and diagonal) substitutions this should lead to wrong solution because refinement isperformed with original matrix. So, refinement stepsare incompatible with forward (and diagonal) solving steps, butDSSprovides full controlof these parameters (MKL_DSS_REFINEMENT_OFF and MKL_DSS_REFINEMENT_ON, MKL_DSS_FORWARD_SOLVE etc.). Of course, PARDISO isoperatedby corresponding iparm values, too.
Probably it will be better to switch off refinement steps in the case of forward and diagonal solving steps, but in this case it will be impossible to switchthem ON in somespecial cases, if any.

With best regards,
Sergey