Add sparskit

sparskit2/BLASSM/README

@@ -0,0 +1,83 @@
+c----------------------------------------------------------------------c
+c                                                                      c
+c                BLASSM  and MATVEC MODULES                            c
+c                                                                      c
+c----------------------------------------------------------------------c
+c                                                                      c 
+c    This directory contains the BLASSM and MATVEC modules of SPARSKIT c
+c                                                                      c
+c----------------------------------------------------------------------c
+c                                                                      c
+c                Current contents                                      c 
+c----------------------------------------------------------------------c
+c                                                                      c 
+c    blassm.f : contains the latest version of the basc linear algerba c
+c               routines for sparse matrices.                          c
+c    								       c
+c    tester.f : is a main program to test the routines and the paths   c
+c                                                                      c
+c    matvec.f : contains the subroutines in the module matvec          c
+c                                                                      c
+c    rmatvec.f: a test program that runs all the routines in matvec    c
+c                                                                      c
+c----------------------------------------------------------------------c
+c                                                                      c
+c    makefile : make file for tester.ex (tests blassm.f) and mvec.ex   c
+c               (tests routines in matvec.f)                           c
+c                                                                      c
+c----------------------------------------------------------------------c
+c----------------------------------------------------------------------c   
+c----------------------------------------------------------------------c
+c              current status of blassm.f                              c
+c                                                                      c
+c----------------------------------------------------------------------c
+c                          S P A R S K I T                             c
+c----------------------------------------------------------------------c
+c        BASIC LINEAR ALGEBRA FOR SPARSE MATRICES. BLASSM MODULE       c
+c----------------------------------------------------------------------c
+c amub   :   computes     C = A*B                                      c
+c aplb   :   computes     C = A+B                                      c
+c aplsb  :   computes     C = A + s B                                  c
+c apmbt  :   Computes     C = A +/- transp(B)                          c
+c aplsbt :   Computes     C = A + s * transp(B)                        c
+c diamua :   Computes     C = Diag * A                                 c
+c amudia :   Computes     C = A* Diag                                  c
+c apldia :   Computes     C = A + Diag.                                c
+c aplsca :   Computes     A:= A + s I    (s = scalar)                  c 
+c----------------------------------------------------------------------c 
+c----------------------------------------------------------------------c
+c                                                                      c
+c              current status of matvec.f                              c
+c                                                                      c
+c----------------------------------------------------------------------c
+c                          S P A R S K I T                             c
+c----------------------------------------------------------------------c
+c          BASIC MATRIX-VECTOR OPERATIONS - MATVEC MODULE              c
+c         Matrix-vector Mulitiplications and Triang. Solves            c
+c----------------------------------------------------------------------c
+c contents:
+c----------                                                            c
+c 1) Matrix-vector products:                                           c
+c---------------------------                                           c
+c amux  : A times a vector. Compressed Sparse Row (CSR) format.        c
+c amuxms: A times a vector. Modified Compress Sparse Row format.       c
+c atmux : Transp(A) times a vector. CSR format.                        c
+c amuxe : A times a vector. Ellpack/Itpack (ELL) format.               c
+c amuxd : A times a vector. Diagonal (DIA) format.                     c
+c amuxj : A times a vector. Jagged Diagonal (JAD) format.              c
+c vbrmv : Sparse matrix-full vector product, in VBR format             c
+c                                                                      c
+c 2) Triangular system solutions:                                      c
+c-------------------------------                                       c
+c lsol  : Unit Lower Triang. solve. Compressed Sparse Row (CSR) format.c
+c ldsol : Lower Triang. solve.  Modified Sparse Row (MSR) format.      c
+c lsolc : Unit Lower Triang. solve. Comp. Sparse Column (CSC) format.  c
+c ldsolc: Lower Triang. solve. Modified Sparse Column (MSC) format.    c
+c ldsoll: Lower Triang. solve with level scheduling. MSR format.       c
+c usol  : Unit Upper Triang. solve. Compressed Sparse Row (CSR) format.c
+c udsol : Upper Triang. solve.  Modified Sparse Row (MSR) format.      c
+c usolc : Unit Upper Triang. solve. Comp. Sparse Column (CSC) format.  c
+c udsolc: Upper Triang. solve.  Modified Sparse Column (MSC) format.   c
+c----------------------------------------------------------------------c
+
+

sparskit2/BLASSM/apl.f

@@ -0,0 +1,91 @@
+      subroutine aplsb (nrow,ncol,a,ja,ia,s,b,jb,ib,c,jc,ic,nzmax,
+     *     iw,ierr)
+      real*8 a(*), b(*), c(*), s
+      integer ja(*),jb(*),jc(*),ia(nrow+1),ib(nrow+1),ic(nrow+1),
+     *     iw(ncol)
+c-----------------------------------------------------------------------
+c performs the matrix sum  C = A+s*B. 
+c-----------------------------------------------------------------------
+c on entry:
+c ---------
+c nrow     = integer. The row dimension of A and B
+c ncol  = integer. The column dimension of A and B.
+c job   = integer. Job indicator. When job = 0, only the structure
+c                  (i.e. the arrays jc, ic) is computed and the
+c                  real values are ignored.
+c
+c a,
+c ja,
+c ia   = Matrix A in compressed sparse row format.
+c 
+c s    = real*8 - coefficient that multiplies B.
+c b, 
+c jb, 
+c ib     =  Matrix B in compressed sparse row format.
+c
+c nzmax	= integer. The  length of the arrays c and jc.
+c         amub will stop if the result matrix C  has a number 
+c         of elements that exceeds exceeds nzmax. See ierr.
+c 
+c on return:
+c----------
+c c, 
+c jc, 
+c ic	= resulting matrix C in compressed sparse row sparse format.
+c	    
+c ierr	= integer. serving as error message. 
+c         ierr = 0 means normal return,
+c         ierr .gt. 0 means that aplsb1 stopped while computing the
+c         i-th row  of C with i=ierr, because the number 
+c         of elements in C exceeds nzmax.
+c
+c work arrays:
+c------------
+c iw	= integer work array of length equal to the number of
+c         columns in A.
+c note: expanded  row implementation. Does not require column indices to
+c       be sorted.
+c-----------------------------------------------------------------------
+      ierr = 0
+      len = 0
+      ic(1) = 1 
+      do 1 j=1, ncol
+         iw(j) = 0
+ 1    continue
+c     
+      do 500 ii=1, nrow
+c     copy row ii  to C
+         do 200 ka=ia(ii), ia(ii+1)-1 
+            len = len+1
+            jcol    = ja(ka)
+            if (len .gt. nzmax) goto 999
+            jc(len) = jcol 
+            c(len)  = a(ka) 
+            iw(jcol)= len
+ 200     continue
+c     
+         do 300 kb=ib(ii),ib(ii+1)-1
+            jcol = jb(kb)
+            jpos = iw(jcol)
+            if (jpos .eq. 0) then
+               len = len+1
+               if (len .gt. nzmax) goto 999
+               jc(len) = jcol
+               c(len)  = s*b(kb)
+               iw(jcol)= len
+            else
+               c(jpos) = c(jpos) + s*b(kb)
+            endif
+ 300     continue
+         do 301 k=ic(ii), len
+         iw(jc(k)) = 0
+ 301     continue
+         ic(ii+1) = len+1
+ 500  continue
+      return
+ 999  ierr = ii
+      return
+c------------end of aplsb1 --------------------------------------------- 
+c-----------------------------------------------------------------------
+      end
+c-----------------------------------------------------------------------