1c96caaccSSatish Balay! 2c96caaccSSatish Balay! Fortran kernel for gemv() BLAS operation. This version supports 3c96caaccSSatish Balay! matrix array stored in single precision but vectors in double 4c96caaccSSatish Balay! 5c96caaccSSatish Balay#include <petsc/finclude/petscsys.h> 60113e719SMartin Diehl 7*0ccf82acSMartin Diehlpure subroutine MSGemv(bs,ncols,A,x,y) 8*0ccf82acSMartin Diehl implicit none (type, external) 9*0ccf82acSMartin Diehl PetscInt, intent(in) :: bs,ncols 10*0ccf82acSMartin Diehl MatScalar, intent(in) :: A(bs,ncols) 11*0ccf82acSMartin Diehl PetscScalar, intent(in) :: x(ncols) 12*0ccf82acSMartin Diehl PetscScalar, intent(out) :: y(bs) 13c96caaccSSatish Balay 14c96caaccSSatish Balay PetscInt i,j 15c96caaccSSatish Balay 160113e719SMartin Diehl do j=1,bs 17c96caaccSSatish Balay y(j) = 0.0d0 180113e719SMartin Diehl end do 19c96caaccSSatish Balay 200113e719SMartin Diehl do i=1,ncols 210113e719SMartin Diehl do j=1,bs 22c96caaccSSatish Balay y(j) = y(j) + A(j,i)*x(i) 230113e719SMartin Diehl end do 240113e719SMartin Diehl end do 250113e719SMartin Diehlend subroutine MSGemv 26c96caaccSSatish Balay 27*0ccf82acSMartin Diehlpure subroutine MSGemvp(bs,ncols,A,x,y) 28*0ccf82acSMartin Diehl implicit none (type, external) 29*0ccf82acSMartin Diehl PetscInt, intent(in) :: bs,ncols 30*0ccf82acSMartin Diehl MatScalar, intent(in) :: A(bs,ncols) 31*0ccf82acSMartin Diehl PetscScalar, intent(in) :: x(ncols) 32*0ccf82acSMartin Diehl PetscScalar, intent(inout) :: y(bs) 33c96caaccSSatish Balay 34c96caaccSSatish Balay PetscInt i, j 35c96caaccSSatish Balay 360113e719SMartin Diehl do i=1,ncols 370113e719SMartin Diehl do j=1,bs 38c96caaccSSatish Balay y(j) = y(j) + A(j,i)*x(i) 390113e719SMartin Diehl end do 400113e719SMartin Diehl end do 410113e719SMartin Diehlend subroutine MSGemvp 42c96caaccSSatish Balay 43*0ccf82acSMartin Diehlpure subroutine MSGemvm(bs,ncols,A,x,y) 44*0ccf82acSMartin Diehl implicit none (type, external) 45*0ccf82acSMartin Diehl PetscInt, intent(in) :: bs,ncols 46*0ccf82acSMartin Diehl MatScalar, intent(in) :: A(bs,ncols) 47*0ccf82acSMartin Diehl PetscScalar, intent(in) :: x(ncols) 48*0ccf82acSMartin Diehl PetscScalar, intent(inout) :: y(bs) 49c96caaccSSatish Balay 50c96caaccSSatish Balay PetscInt i, j 51c96caaccSSatish Balay 520113e719SMartin Diehl do i=1,ncols 530113e719SMartin Diehl do j=1,bs 54c96caaccSSatish Balay y(j) = y(j) - A(j,i)*x(i) 550113e719SMartin Diehl end do 560113e719SMartin Diehl end do 570113e719SMartin Diehlend subroutine MSGemvm 58c96caaccSSatish Balay 59*0ccf82acSMartin Diehlpure subroutine MSGemvt(bs,ncols,A,x,y) 60*0ccf82acSMartin Diehl implicit none (type, external) 61*0ccf82acSMartin Diehl PetscInt, intent(in) :: bs,ncols 62*0ccf82acSMartin Diehl MatScalar, intent(in) :: A(bs,ncols) 63*0ccf82acSMartin Diehl PetscScalar, intent(in) :: x(bs) 64*0ccf82acSMartin Diehl PetscScalar, intent(inout) :: y(ncols) 65c96caaccSSatish Balay 66c96caaccSSatish Balay PetscInt i,j 67c96caaccSSatish Balay PetscScalar sum 68*0ccf82acSMartin Diehl 690113e719SMartin Diehl do i=1,ncols 70c96caaccSSatish Balay sum = y(i) 710113e719SMartin Diehl do j=1,bs 72c96caaccSSatish Balay sum = sum + A(j,i)*x(j) 730113e719SMartin Diehl end do 74c96caaccSSatish Balay y(i) = sum 750113e719SMartin Diehl end do 760113e719SMartin Diehlend subroutine MSGemvt 77c96caaccSSatish Balay 78*0ccf82acSMartin Diehlpure subroutine MSGemm(bs,A,B,C) 79*0ccf82acSMartin Diehl implicit none (type, external) 80*0ccf82acSMartin Diehl PetscInt, intent(in) :: bs 81*0ccf82acSMartin Diehl MatScalar, intent(in) :: B(bs,bs),C(bs,bs) 82*0ccf82acSMartin Diehl MatScalar, intent(inout) :: A(bs,bs) 83*0ccf82acSMartin Diehl 84c96caaccSSatish Balay PetscScalar sum 85c96caaccSSatish Balay PetscInt i,j,k 86c96caaccSSatish Balay 870113e719SMartin Diehl do i=1,bs 880113e719SMartin Diehl do j=1,bs 89c96caaccSSatish Balay sum = A(i,j) 900113e719SMartin Diehl do k=1,bs 91c96caaccSSatish Balay sum = sum - B(i,k)*C(k,j) 920113e719SMartin Diehl end do 93c96caaccSSatish Balay A(i,j) = sum 940113e719SMartin Diehl end do 950113e719SMartin Diehl end do 960113e719SMartin Diehlend subroutine MSGemm 97c96caaccSSatish Balay 98*0ccf82acSMartin Diehlpure subroutine MSGemmi(bs,A,C,B) 99*0ccf82acSMartin Diehl implicit none (type, external) 100*0ccf82acSMartin Diehl PetscInt, intent(in) :: bs 101*0ccf82acSMartin Diehl MatScalar, intent(in) :: B(bs,bs),C(bs,bs) 102*0ccf82acSMartin Diehl MatScalar, intent(out) :: A(bs,bs) 103c96caaccSSatish Balay 104*0ccf82acSMartin Diehl PetscScalar sum 105c96caaccSSatish Balay PetscInt i,j,k 106c96caaccSSatish Balay 1070113e719SMartin Diehl do i=1,bs 1080113e719SMartin Diehl do j=1,bs 109c96caaccSSatish Balay sum = 0.0d0 1100113e719SMartin Diehl do k=1,bs 111c96caaccSSatish Balay sum = sum + B(i,k)*C(k,j) 1120113e719SMartin Diehl end do 113c96caaccSSatish Balay A(i,j) = sum 1140113e719SMartin Diehl end do 1150113e719SMartin Diehl end do 1160113e719SMartin Diehlend subroutine MSGemmi 117