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