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:       PetscBool   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 MatCreate(PETSC_COMM_WORLD,A,ierr)
 28:       call MatSetType(A, MATMPIAIJ,ierr)
 29:       call MatLoad(A,v,ierr)
 30:       CHKERRQ(ierr)
 31:       call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)

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

 52: !
 53: !   Print off-diagonal portion of matrix
 54: !
 55:       call MatGetRowIJ(Ao,one,zero,zero,n,ia,iia,ja,jja,done,ierr)
 56:       call MatGetArray(Ao,aa,aaa,ierr)
 57:       do 30, i=1,n
 58:         write(7+rank,*) 'row ',i+rstart,' number nonzeros ',               &
 59:      &                  ia(iia+i+1)-ia(iia+i)
 60:         do 40, j=ia(iia+i),ia(iia+i+1)-1
 61:           write(7+rank,*)'  ',j,icol(iicol+ja(jja+j))+1,aa(aaa+j)
 62:  40     continue
 63:  30   continue
 64:       call MatRestoreRowIJ(Ao,one,zero,zero,n,ia,iia,ja,jja,done,ierr)
 65:       call PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr)

 67:       call MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr)
 68:       call MatDestroy(A,ierr)
 69:       call PetscViewerDestroy(v,ierr)
 70:       call PetscFinalize(ierr)
 71:       end