Actual source code: sgemv.F90

  1: !
  2: !    Fortran kernel for gemv() BLAS operation. This version supports
  3: !  matrix array stored in single precision but vectors in double
  4: !
  5: #include <petsc/finclude/petscsys.h>
  6: !
  7:       subroutine MSGemv(bs,ncols,A,x,y)
  8:       implicit none
  9:       PetscInt          bs,ncols
 10:       MatScalar        A(bs,ncols)
 11:       PetscScalar      x(ncols),y(bs)

 13:       PetscInt         i,j

 15:       do 5, j=1,bs
 16:         y(j) = 0.0d0
 17:  5    continue

 19:       do 10, i=1,ncols
 20:         do 20, j=1,bs
 21:           y(j) = y(j) + A(j,i)*x(i)
 22:  20     continue
 23:  10   continue

 25:       end

 27:       subroutine MSGemvp(bs,ncols,A,x,y)
 28:       implicit none
 29:       PetscInt          bs,ncols
 30:       MatScalar        A(bs,ncols)
 31:       PetscScalar      x(ncols),y(bs)

 33:       PetscInt         i, j

 35:       do 10, i=1,ncols
 36:         do 20, j=1,bs
 37:           y(j) = y(j) + A(j,i)*x(i)
 38:  20     continue
 39:  10   continue

 41:       end

 43:       subroutine MSGemvm(bs,ncols,A,x,y)
 44:       implicit none
 45:       PetscInt          bs,ncols
 46:       MatScalar        A(bs,ncols)
 47:       PetscScalar      x(ncols),y(bs)

 49:       PetscInt         i, j

 51:       do 10, i=1,ncols
 52:         do 20, j=1,bs
 53:           y(j) = y(j) - A(j,i)*x(i)
 54:  20     continue
 55:  10   continue

 57:       end

 59:       subroutine MSGemvt(bs,ncols,A,x,y)
 60:       implicit none
 61:       PetscInt          bs,ncols
 62:       MatScalar        A(bs,ncols)
 63:       PetscScalar      x(bs),y(ncols)

 65:       PetscInt          i,j
 66:       PetscScalar      sum
 67:       do 10, i=1,ncols
 68:         sum = y(i)
 69:         do 20, j=1,bs
 70:           sum = sum + A(j,i)*x(j)
 71:  20     continue
 72:         y(i) = sum
 73:  10   continue

 75:       end

 77:       subroutine MSGemm(bs,A,B,C)
 78:       implicit none
 79:       PetscInt    bs
 80:       MatScalar   A(bs,bs),B(bs,bs),C(bs,bs)
 81:       PetscScalar sum
 82:       PetscInt    i,j,k

 84:       do 10, i=1,bs
 85:         do 20, j=1,bs
 86:           sum = A(i,j)
 87:           do 30, k=1,bs
 88:             sum = sum - B(i,k)*C(k,j)
 89:  30       continue
 90:           A(i,j) = sum
 91:  20     continue
 92:  10   continue

 94:       end

 96:       subroutine MSGemmi(bs,A,C,B)
 97:       implicit none
 98:       PetscInt    bs
 99:       MatScalar   A(bs,bs),B(bs,bs),C(bs,bs)
100:       PetscScalar sum

102:       PetscInt    i,j,k

104:       do 10, i=1,bs
105:         do 20, j=1,bs
106:           sum = 0.0d0
107:           do 30, k=1,bs
108:             sum = sum + B(i,k)*C(k,j)
109:  30       continue
110:           A(i,j) = sum
111:  20     continue
112:  10   continue

114:       end