Actual source code: sgemv.F

  1: !
  2: !  "$Id: sgemv.F,v 1.4 2001/08/07 03:05:24 balay Exp $";
  3: !
  4: !    Fortran kernel for gemv() BLAS operation. This version supports
  5: !  matrix array stored in single precision but vectors in double
  6: !
 7:  #include include/finclude/petscdef.h
  8: !
  9:       subroutine MSGemv(bs,ncols,A,x,y)
 10:       implicit none
 11:       integer          bs,ncols
 12:       MatScalar        a(bs,ncols)
 13:       PetscScalar      x(ncols),y(bs)

 15:       integer          i,j

 17:       do 5, j=1,bs
 18:         y(j) = 0.0d0
 19:  5    continue

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

 27:       return
 28:       end

 30:       subroutine MSGemvp(bs,ncols,A,x,y)
 31:       implicit none
 32:       integer          bs,ncols
 33:       MatScalar        a(bs,ncols)
 34:       PetscScalar      x(ncols),y(bs)

 36:       integer          i,j

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

 44:       return
 45:       end

 47:       subroutine MSGemvm(bs,ncols,A,x,y)
 48:       implicit none
 49:       integer          bs,ncols
 50:       MatScalar        a(bs,ncols)
 51:       PetscScalar      x(ncols),y(bs)

 53:       integer          i,j

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

 61:       return
 62:       end

 64:       subroutine MSGemvt(bs,ncols,A,x,y)
 65:       implicit none
 66:       integer          bs,ncols
 67:       MatScalar        a(bs,ncols)
 68:       PetscScalar      x(bs),y(ncols)

 70:       integer          i,j
 71:       PetscScalar      sum

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

 81:       return
 82:       end

 84:       subroutine MSGemm(bs,A,B,C)
 85:       implicit none
 86:       integer          bs
 87:       MatScalar        A(bs,bs),B(bs,bs),C(bs,bs)
 88:       PetscScalar      sum

 90:       integer          i,j,k

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

102:       return
103:       end

105:       subroutine MSGemmi(bs,A,C,B)
106:       implicit none
107:       integer          bs
108:       MatScalar        A(bs,bs),B(bs,bs),C(bs,bs)
109:       PetscScalar      sum

111:       integer          i,j,k

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

123:       return
124:       end