Actual source code: ex79f.F
1: !
2: ! This program demonstrates use of MatGetRowIJ() from Fortran
3: !
4: program main
5: implicit none
6: #include finclude/petscsys.h
7: #include finclude/petscmat.h
8: #include finclude/petscis.h
9: #include finclude/petscviewer.h
11: Mat A,Ad,Ao
12: PetscErrorCode ierr
13: PetscMPIInt rank
14: PetscViewer v
15: PetscInt i,j,ia(1),ja(1)
16: PetscInt n,icol(1),rstart
17: PetscInt zero,one,rend
18: PetscTruth done
19: PetscOffset iia,jja,aaa,iicol
20: PetscScalar aa(1)
22: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
23: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
25: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,'../matbinary.ex', &
26: & FILE_MODE_READ,v,ierr)
27: call MatLoad(v,MATMPIAIJ,A,ierr)
28: CHKERRQ(ierr)
29: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
31: call MatMPIAIJGetSeqAIJ(A,Ad,Ao,icol,iicol,ierr)
32: call MatGetOwnershipRange(A,rstart,rend,ierr)
33: !
34: ! Print diagonal portion of matrix
35: !
36: call PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr)
37: zero = 0
38: one = 1
39: call MatGetRowIJ(Ad,one,zero,zero,n,ia,iia,ja,jja,done,ierr)
40: call MatGetArray(Ad,aa,aaa,ierr)
41: do 10, i=1,n
42: write(7+rank,*) 'row ',i+rstart,' number nonzeros ', &
43: & ia(iia+i+1)-ia(iia+i)
44: do 20, j=ia(iia+i),ia(iia+i+1)-1
45: write(7+rank,*)' ',j,ja(jja+j)+rstart,aa(aaa+j)
46: 20 continue
47: 10 continue
48: call MatRestoreRowIJ(Ad,one,zero,zero,n,ia,iia,ja,jja,done,ierr)
50: !
51: ! Print off-diagonal portion of matrix
52: !
53: call MatGetRowIJ(Ao,one,zero,zero,n,ia,iia,ja,jja,done,ierr)
54: call MatGetArray(Ao,aa,aaa,ierr)
55: do 30, i=1,n
56: write(7+rank,*) 'row ',i+rstart,' number nonzeros ', &
57: & ia(iia+i+1)-ia(iia+i)
58: do 40, j=ia(iia+i),ia(iia+i+1)-1
59: write(7+rank,*)' ',j,icol(iicol+ja(jja+j))+1,aa(aaa+j)
60: 40 continue
61: 30 continue
62: call MatRestoreRowIJ(Ao,one,zero,zero,n,ia,iia,ja,jja,done,ierr)
63: call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr)
65: call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
66: call MatDestroy(A,ierr)
67: call PetscViewerDestroy(v,ierr)
68: call PetscFinalize(ierr)
69: end