topic Here is my codes. in Intel® oneAPI Math Kernel Library

dss_solver gives non-stable solve

Tianyu_D_ — Mon, 07 Apr 2014 15:44:37 GMT

Hi there,

I am using the mkl_dss.f77 to solve large sparse non-linear equations by iteratively calling the dss_solver. I just spotted that although I set the same initial values, I could not obtain a same result for every time. And I found this is due to the dss_solver gives results with slightly difference every time and the difference adds up in the iteration process which leads to a non-stable solve for the large sparse non-linear equations.

How can I make the mkl_dss solver to have repetitive accuracy? Anyone could give me some advice on this subject?

Hi,

Alexander_K_Intel2 — Mon, 07 Apr 2014 16:36:22 GMT

Hi,

It is better to use Pardiso interface for this aim. In Pardiso interface you can use CNR mode of MKL or set iparm[33] to maximum number of threads - both approach will provide stable results from run to run.

Thanks,

Alex

Hi Alex,

Tianyu_D_ — Mon, 07 Apr 2014 16:55:57 GMT

Hi Alex,

Thank you for your reply.

The reason I did not use Pardiso is that DSS provides the function dss_reorder which could help to save time since I only need to reorder once in the iteration process.

So if I want to keep using DSS, any suggestion for me to make its output stable?

Regards,

Tianyu

and also matrix checker

Gennady_F_Intel — Mon, 07 Apr 2014 16:59:27 GMT

and also matrix checker option available with Pardiso API too. see ipatr(27). This functionality is not available for DSS API.

Hi Tianyu,

Alexander_K_Intel2 — Mon, 07 Apr 2014 17:06:47 GMT

Hi Tianyu,

Almost all that can be done via DSS interface can be implemented via PARDISO. Reordering step DSS_reorder correspondent to pardiso call with ido=11. But it is not quite good approach to call reordering only once for all matrix except symmetric-positive define, because good reordering step for all other matrix significantly depend on matrix value. Can you provide some pseudo-code of you algorithm that we can change to show how it can be done via PARDISO interface?

Thanks,

Alex

Here is my codes.

Tianyu_D_ — Mon, 07 Apr 2014 17:20:34 GMT

Here is my codes.

Thank you so much.

 SUBROUTINE NFDSS_SOLVE(XDIFF,IU,JU,UM,nRow0,ITER,Idum)
      INCLUDE 'mkl_dss.f77'
      !INCLUDE 'mkl_pardiso.f77'

      INTEGER IU(*), JU(*), Idum(*), nRow0, nRow, ITER
      REAL  UM(*), XDIFF(*)
      
      INTEGER tIU(:), tJU(:)
      DOUBLE PRECISION tUM(:), tXDIFF(:), tX(:)
      ALLOCATABLE tIU, tJU, tUM, tXDIFF, tX
      
      INTEGER nNonZeros, nNonZeros0
      
      INTEGER*8 handle
      INTEGER error

      nNonZeros0=IU(nRow0+1)-1;
      
      nRow=nRow0-2;
      nNonZeros=nNonZeros0-2;
      
      ALLOCATE(tIU(nRow+1))
      ALLOCATE(tXDIFF(nRow))
      ALLOCATE(tX(nRow))
      ALLOCATE(tJU(nNonZeros))
      ALLOCATE(tUM(nNonZeros))
      
      tIU(1:nRow+1)=IU(3:nRow+3)-2
      tJU(1:nNonZeros)=JU(3:nNonZeros+2)-2
      tUM(1:nNonZeros)=UM(3:nNonZeros+2)
      tXDIFF(1:nRow)=XDIFF(3:nRow+3)

    
C---------------------------------------------------------------------------
C Initialize the solver
C---------------------------------------------------------------------------
      error = dss_create(handle, MKL_DSS_DEFAULTS+MKL_DSS_OOC_VARIABLE)
      IF (error .NE. MKL_DSS_SUCCESS ) GOTO 999
C---------------------------------------------------------------------------

C---------------------------------------------------------------------------
C Define the non-zero structure of the matrix
C---------------------------------------------------------------------------
      error = dss_define_structure( handle, MKL_DSS_NON_SYMMETRIC,
     & tIU, nRow, nRow, tJU, nNonZeros )
      IF (error .NE. MKL_DSS_SUCCESS ) GOTO 999
C---------------------------------------------------------------------------

C---------------------------------------------------------------------------
C Reorder the matrix
C---------------------------------------------------------------------------

      IF (ITER.EQ.1) THEN
          error = dss_reorder( handle, MKL_DSS_GET_ORDER, Idum)
          IF (error .NE. MKL_DSS_SUCCESS ) GOTO 999
      ELSE
          error = dss_reorder( handle, MKL_DSS_MY_ORDER, Idum)
          IF (error .NE. MKL_DSS_SUCCESS ) GOTO 999
      ENDIF

C---------------------------------------------------------------------------
C Factor the matrix
C---------------------------------------------------------------------------
      error = dss_factor_real( handle, MKL_DSS_DEFAULTS, tUM)
      IF (error .NE. MKL_DSS_SUCCESS ) GOTO 999
      
C---------------------------------------------------------------------------
C Get the solution vector
C---------------------------------------------------------------------------
      error = dss_solve_real( handle, MKL_DSS_DEFAULTS,
     & tXDIFF, 1, tX)
      IF (error .NE. MKL_DSS_SUCCESS ) GOTO 999
C---------------------------------------------------------------------------


      
C---------------------------------------------------------------------------
C Deallocate solver storage
C---------------------------------------------------------------------------
      error = dss_delete( handle, MKL_DSS_DEFAULTS )
      IF (error .NE. MKL_DSS_SUCCESS ) GOTO 999
C---------------------------------------------------------------------------
 
      GOTO 1999
  999 WRITE(*,*) "Solver returned error code ", error
      
 1999 XDIFF(3:nRow+2)= tX(1:nRow)
      
      
      DEALLOCATE(tIU,tJU,tUM,tXDIFF,tX)
      RETURN
      END

Hi Alex,

Tianyu_D_ — Mon, 07 Apr 2014 17:27:03 GMT

Hi Alex,

The attached is my code.