Actual source code: sgemv.F

  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 include/finclude/petscdef.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:       return
 26:       end

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

 34:       PetscInt          i,j

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

 42:       return
 43:       end

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

 51:       PetscInt          i,j

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

 59:       return
 60:       end

 62:       subroutine MSGemvt(bs,ncols,A,x,y)
 63:       implicit none
 64:       PetscInt          bs,ncols
 65:       MatScalar        A(bs,ncols)
 66:       PetscScalar      x(bs),y(ncols)

 68:       PetscInt          i,j
 69:       PetscScalar      sum

 71:       do 10, i=1,ncols
 72:         sum = y(i)
 73:         do 20, j=1,bs
 74:           sum = sum + A(j,i)*x(j)
 75:  20     continue
 76:         y(i) = sum
 77:  10   continue

 79:       return
 80:       end

 82:       subroutine MSGemm(bs,A,B,C)
 83:       implicit none
 84:       PetscInt   bs
 85:       MatScalar  A(bs,bs),B(bs,bs),C(bs,bs)
 86:       PetscScalar sum

 88:       PetscInt   i,j,k

 90:       do 10, i=1,bs
 91:         do 20, j=1,bs
 92:           sum = A(i,j)
 93:           do 30, k=1,bs
 94:             sum = sum - B(i,k)*C(k,j)
 95:  30       continue
 96:           A(i,j) = sum
 97:  20     continue
 98:  10   continue

100:       return
101:       end

103:       subroutine MSGemmi(bs,A,C,B)
104:       implicit none
105:       PetscInt    bs
106:       MatScalar   A(bs,bs),B(bs,bs),C(bs,bs)
107:       PetscScalar sum

109:       PetscInt    i,j,k

111:       do 10, i=1,bs
112:         do 20, j=1,bs
113:           sum = 0.0d0
114:           do 30, k=1,bs
115:             sum = sum + B(i,k)*C(k,j)
116:  30       continue
117:           A(i,j) = sum
118:  20     continue
119:  10   continue

121:       return
122:       end