CSR MATRIX + ARPACK

[vickysharma0812]

Add support to find out the Eigenvectors and Eigenvalues of sparse matrix.
To this end use ARPACK library.
The library is here
https://www.caam.rice.edu/software/ARPACK/

• Symmetric matrix smallest eigenvalue
• Symmetric matrix largest eigenvalue
• Symmetric matrix smallest eigenvalue and eigenvector
• Symmetric matrix largest eigenvalue and eigenvector
• Generalized eigenvalue problem
• Eigenvalue and eigenvectors of generalized matrix

About ARPACK:

ARPACK is a collection of Fortran77 subroutines designed to solve large scale eigenvalue problems.

• The package is designed to compute a few eigenvalues and corresponding eigenvectors of a general n by n matrix A.
• ARPACK software is capable of solving large scale symmetric, nonsymmetric, and generalized eigenproblems from significant application areas.
• For many standard problems, a matrix factorization is not required. Only the action of the matrix on a vector is needed.

You can download the software from here : https://www.caam.rice.edu/software/ARPACK/download.html

Regards
Vikas

[vickysharma0812]

References

See this

• FIDASIM

Add more details

[vickysharma0812]

References

See this

• utilities@FIDASIM

[vickysharma0812]

The work on this issue started on 2022-12-08.

[vickysharma0812]

Progress made so far:

Import modules

PROGRAM main
USE easifemBase
Implicit none

Understanding the arguments

Arguments to SAUPD

Argument	Type	Intent	Value
IDO	Int	INOUT	0
BMAT	Char(1)	IN	"I”
N	Int	IN	SIZE(A,1)
WHICH	Char(2)	IN	LA
NEV	Int	IN	1
TOL	Real32	IN	1.0E-10
RESID	Real32(N)	INOUT	NA, as info .EQ. 0
NCV	Int	IN	3*N
V	Real32(N, NCV)	OUT
LDV	Int	IN	SIZE(V,1)
IPARAM	Int (11)	INOUT	See below
IPNTR	Int (11)	OUT	NA
WORKD	Real32(3*N)	INOUT	NA
WORKL	Real32(LWORKL)	OUT	NA
LWORKL	Int	IN	at least NCV**2 + 8*NCV
INFO	Int	INOUT	0

Arguments to IPARAM:

		Options	Selected value	Intent
IPARAM (1)	ISHIFT	1,2	1	IN
IPARAM (2)	LEVEC	NA	NA	IN: Deprecated
IPARAM (3)	MAXITER		3*N	IN
IPARAM (4)	NB	1	1	IN
IPARAM (5)	NCONV			OUT
IPARAM (6)	IUPD	NA	NA	IN: Deprecated
IPARAM (7)	MODE	1,2,3,4,5	1	IN
IPARAM (8)	NP		NA	IN:
IPARAM (9)	NUMOP			OUT
IPARAM (10)	NUMOPB			OUT
IPARAM (11)	NUMREO			OUT

Declare variables

REAL( DFP ) :: maxEV
REAL( DFP ), ALLOCATABLE :: mat(:,:)
INTEGER( I4B ) :: ncv

Getting the algebraic largest eigenvalue of a diagonal matrix.

  mat = zeros(100,100, 1.0_DFP)
  call SetDiag(mat=mat, d=arange(1, SIZE(mat,1)), diagNo=0)
  maxEV = SymLargestEigenVal(mat=mat)
  CALL Display(maxEV, "maxEV=")

Getting the absolute largest eigenvalue of a diagonal matrix. In this case we
provide extra argument which="LM".

:::note Default value of which is "LA" which stands for largest ALGEBRAIC eigenvalue.
:::

  call SetDiag(mat=mat, d=arange(1, SIZE(mat,1)), diagNo=0)
  mat(SIZE(mat,1), SIZE(mat,2) ) = -1000
  maxEV = SymLargestEigenVal(mat=mat, which="LM" )
  CALL Display(maxEV, "max absolute EV=")

:::caution When which="LA", the returned eigenvalue can be positive.
:::

END PROGRAM main

[vickysharma0812]

Eigenvalue of diffusion matrix with zero dirichlet boundary condition.

Import modules

program main
use easifembase
implicit none

Declare variables

  real( dfp ) :: maxev
  real( dfp ), allocatable :: mat(:,:)
  integer( i4b ) :: ncv
  integer(I4B) :: nx, ny, n, ii
  integer(I4B), allocatable :: dbcnptrs(:)
  real(DFP) :: inv_dx2

  nx = 10; ny=nx; n = nx*ny
  inv_dx2 = 1.0 !REAL(nx-1, kind=DFP) * REAL(nx-1, kind=DFP)

  mat = zeros(n,n, 1.0_DFP)
  call SetDiag(mat=mat, d=inv_dx2 * [4.0], diagNo=0)
  call SetDiag(mat=mat, d=inv_dx2 * [-1.0], diagNo=1)
  call SetDiag(mat=mat, d=inv_dx2 * [-1.0], diagNo=5)
  call GetSym(mat=mat, from="U")

  call display( MdEncode(mat(1:5, 1:5)), "mat(1:5, 1:5) = " )

  dbcnptrs = (arange(2, nx-1, 1) .APPEND. arange(n-nx+2, n-1, 1)) &
    & .APPEND. &
    & (arange(1, nx*ny, nx ) .APPEND. arange(nx, nx*(ny-1)+1, nx))

  mat(:, dbcnptrs) = 0.0_DFP
  mat(dbcnptrs, :) = 0.0_DFP
  do concurrent(ii=1:size(dbcnptrs))
    mat(dbcnptrs(ii), dbcnptrs(ii)) = 1.0_DFP
  end do

  call display( MdEncode(mat(1:5, 1:5)), "mat(1:5, 1:5) = " )

  maxEV = SymLargestEigenVal(mat=mat, which="LM")
  call Display(maxEV, "maxEV=")

results

mat(1:5, 1:5) =

4	-1	0	0	0
-1	4	-1	0	0
0	-1	4	-1	0
0	0	-1	4	-1
0	0	0	-1	4

mat(1:5, 1:5) =

1	0	0	0	0
0	1	0	0	0
0	0	1	0	0
0	0	0	1	0
0	0	0	0	1

maxEV=7.48797

end program main

[vickysharma0812]

This example shows the usage of SymLargestEigenValue.

In this example we will calculate many largest eigenvalues.

Import modules

PROGRAM main
USE easifemBase
Implicit none

Declare variables

REAL( DFP ), ALLOCATABLE :: maxEV(:), mat(:,:)
INTEGER( I4B ) :: ncv

Getting the first five algebraic largest eigenvalue of a diagonal matrix.

  mat = zeros(100,100, 1.0_DFP)
  call SetDiag(mat=mat, d=arange(1, SIZE(mat,1)), diagNo=0)
  maxEV = SymLargestEigenVal(mat=mat, nev=5)
  CALL Display(maxEV, "maxEV=")

 maxEV=
-------
 96.000
 97.000
 98.000
 99.000
100.000

Getting the first five absolute largest eigenvalue of a diagonal matrix. In this case we
provide extra argument which="LM".

:::note Default value of which is "LA" which stands for largest absolute eigenvalue.
:::

  call SetDiag(mat=mat, d=arange(1, SIZE(mat,1)), diagNo=0)
  mat(SIZE(mat,1), SIZE(mat,2) ) = -1000
  maxEV = SymLargestEigenVal(mat=mat, nev=5, which="LM" )
  CALL Display(maxEV, "max absolute EV=")

max absolute EV=
----------------
       96.00
       97.00
       98.00
       99.00
    -1000.00

:::caution When which="LA", the returned eigenvalue can be positive.
:::

END PROGRAM main