Actual source code: fmult.F
1: !
2: !
3: ! Fortran kernel for sparse matrix-vector product in the AIJ matrix format
4: !
5: #include include/finclude/petscdef.h
6: !
7: subroutine FortranMultTransposeAddAIJ(n,x,ii,jj,a,y)
8: implicit none
9: PetscScalar x(0:*),a(0:*),y(0:*),alpha
10: PetscInt n,ii(*),jj(0:*)
12: PetscInt i,j,jstart,jend
14: jend = ii(1)
15: do 10,i=1,n
16: jstart = jend
17: jend = ii(i+1)
18: alpha = x(i-1)
19: do 20 j=jstart,jend-1
20: y(jj(j)) = y(jj(j)) + alpha*a(j)
21: 20 continue
22: 10 continue
24: return
25: end
27: subroutine FortranMultAIJ(n,x,ii,jj,a,y)
28: implicit none
29: PetscScalar x(0:*),a(0:*),y(*)
30: PetscInt n,ii(*),jj(0:*)
32: PetscInt i,j,jstart,jend
33: PetscScalar sum
35: jend = ii(1)
36: do 10,i=1,n
37: jstart = jend
38: jend = ii(i+1)
39: sum = 0.d0
40: do 20 j=jstart,jend-1
41: sum = sum + a(j)*x(jj(j))
42: 20 continue
43: y(i) = sum
44: 10 continue
46: return
47: end