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>

  7: pure subroutine MSGemv(bs, ncols, A, x, y)
  8:   use, intrinsic :: ISO_C_binding
  9:   implicit none(type, external)
 10:   PetscInt, intent(in) :: bs, ncols
 11:   MatScalar, intent(in) :: A(bs, ncols)
 12:   PetscScalar, intent(in) :: x(ncols)
 13:   PetscScalar, intent(out) :: y(bs)

 15:   PetscInt :: i

 17:   y(1:bs) = 0.0
 18:   do i = 1, ncols
 19:     y(1:bs) = y(1:bs) + A(1:bs, i)*x(i)
 20:   end do
 21: end subroutine MSGemv

 23: pure subroutine MSGemvp(bs, ncols, A, x, y)
 24:   use, intrinsic :: ISO_C_binding
 25:   implicit none(type, external)
 26:   PetscInt, intent(in) :: bs, ncols
 27:   MatScalar, intent(in) :: A(bs, ncols)
 28:   PetscScalar, intent(in) :: x(ncols)
 29:   PetscScalar, intent(inout) :: y(bs)

 31:   PetscInt :: i

 33:   do i = 1, ncols
 34:     y(1:bs) = y(1:bs) + A(1:bs, i)*x(i)
 35:   end do
 36: end subroutine MSGemvp

 38: pure subroutine MSGemvm(bs, ncols, A, x, y)
 39:   use, intrinsic :: ISO_C_binding
 40:   implicit none(type, external)
 41:   PetscInt, intent(in) :: bs, ncols
 42:   MatScalar, intent(in) :: A(bs, ncols)
 43:   PetscScalar, intent(in) :: x(ncols)
 44:   PetscScalar, intent(inout) :: y(bs)

 46:   PetscInt :: i

 48:   do i = 1, ncols
 49:     y(1:bs) = y(1:bs) - A(1:bs, i)*x(i)
 50:   end do
 51: end subroutine MSGemvm

 53: pure subroutine MSGemvt(bs, ncols, A, x, y)
 54:   use, intrinsic :: ISO_C_binding
 55:   implicit none(type, external)
 56:   PetscInt, intent(in) :: bs, ncols
 57:   MatScalar, intent(in) :: A(bs, ncols)
 58:   PetscScalar, intent(in) :: x(bs)
 59:   PetscScalar, intent(inout) :: y(ncols)

 61:   PetscInt :: i

 63:   do i = 1, ncols
 64:     y(i) = y(i) + sum(A(1:bs, i)*x(1:bs))
 65:   end do
 66: end subroutine MSGemvt

 68: pure subroutine MSGemm(bs, A, B, C)
 69:   use, intrinsic :: ISO_C_binding
 70:   implicit none(type, external)
 71:   PetscInt, intent(in) :: bs
 72:   MatScalar, intent(in) :: B(bs, bs), C(bs, bs)
 73:   MatScalar, intent(inout) :: A(bs, bs)

 75:   PetscInt :: i, j

 77:   do i = 1, bs
 78:     do j = 1, bs
 79:       A(i, j) = A(i, j) - sum(B(i, 1:bs)*C(1:bs, j))
 80:     end do
 81:   end do
 82: end subroutine MSGemm

 84: pure subroutine MSGemmi(bs, A, C, B)
 85:   use, intrinsic :: ISO_C_binding
 86:   implicit none(type, external)
 87:   PetscInt, intent(in) :: bs
 88:   MatScalar, intent(in) :: B(bs, bs), C(bs, bs)
 89:   MatScalar, intent(out) :: A(bs, bs)

 91:   PetscInt :: i, j

 93:   do i = 1, bs
 94:     do j = 1, bs
 95:       A(i, j) = sum(B(i, 1:bs)*C(1:bs, j))
 96:     end do
 97:   end do
 98: end subroutine MSGemmi