Actual source code: ex79f.F

  2: !      "$Id: ex79f.F,v 1.3 2001/08/07 03:03:07 balay Exp $";
  3: !
  4: !   This program demonstrates use of MatGetRowIJ() from Fortran
  5: !
  6:       program main
  7:       implicit none
 8:  #include include/finclude/petsc.h
 9:  #include include/finclude/petscmat.h
 10:  #include include/finclude/petscis.h
 11:  #include include/finclude/petscviewer.h

 13:       Mat         A,Ad,Ao
 14:       integer     ierr,rank
 15:       PetscViewer v
 16:       integer     i,j,ia(1),ja(1),n,icol(1),rstart,rend
 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:       call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)

 29:       call MatMPIAIJGetSeqAIJ(A,Ad,Ao,icol,iicol,ierr)
 30:       call MatGetOwnershipRange(A,rstart,rend,ierr)
 31: !
 32: !   Print diagonal portion of matrix
 33: !
 34:       call PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr)
 35:       call MatGetRowIJ(Ad,1,0,n,ia,iia,ja,jja,done,ierr)
 36:       call MatGetArray(Ad,aa,aaa,ierr)
 37:       do 10, i=1,n
 38:         write(7+rank,*) 'row ',i+rstart,' number nonzeros ',                &
 39:      &                   ia(iia+i+1)-ia(iia+i)
 40:         do 20, j=ia(iia+i),ia(iia+i+1)-1
 41:           write(7+rank,*)'  ',j,ja(jja+j)+rstart,aa(aaa+j)
 42:  20     continue
 43:  10   continue
 44:       call MatRestoreRowIJ(Ad,1,0,n,ia,iia,ja,jja,done,ierr)

 46: !
 47: !   Print off-diagonal portion of matrix
 48: !
 49:       call MatGetRowIJ(Ao,1,0,n,ia,iia,ja,jja,done,ierr)
 50:       call MatGetArray(Ao,aa,aaa,ierr)
 51:       do 30, i=1,n
 52:         write(7+rank,*) 'row ',i+rstart,' number nonzeros ',               &
 53:      &                  ia(iia+i+1)-ia(iia+i)
 54:         do 40, j=ia(iia+i),ia(iia+i+1)-1
 55:           write(7+rank,*)'  ',j,icol(iicol+ja(jja+j))+1,aa(aaa+j)
 56:  40     continue
 57:  30   continue
 58:       call MatRestoreRowIJ(Ao,1,0,n,ia,iia,ja,jja,done,ierr)
 59:       call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr)

 61:       call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
 62:       call MatDestroy(A,ierr)
 63:       call PetscViewerDestroy(v,ierr)
 64:       call PetscFinalize(ierr)
 65:       end