Actual source code: ex79f.F90

  1: !
  2: !   This program demonstrates use of MatGetRowIJ() from Fortran
  3: !
  4:       program main

  6: #include <petsc/finclude/petscmat.h>
  7:       use petscmat
  8:       implicit none

 10:       Mat         A,Ad,Ao
 11:       PetscErrorCode ierr
 12:       PetscMPIInt rank
 13:       PetscViewer v
 14:       PetscInt i,j
 15:       PetscInt n,rstart
 16:       PetscInt zero,one,rend
 17:       PetscBool   done,bb
 18:       PetscScalar,pointer :: aa(:)
 19:       PetscInt,pointer :: ia(:),ja(:),icol(:)

 21:       PetscCallA(PetscInitialize(ierr))
 22:       PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr))

 24:       PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD,'${PETSC_DIR}/share/petsc/datafiles/matrices/' // 'ns-real-int32-float64',FILE_MODE_READ,v,ierr))
 25:       PetscCallA(MatCreate(PETSC_COMM_WORLD,A,ierr))
 26:       PetscCallA(MatSetType(A, MATMPIAIJ,ierr))
 27:       PetscCallA(MatLoad(A,v,ierr))
 28:       PetscCallA(MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr))

 30:       PetscCallA(MatMPIAIJGetSeqAIJF90(A,Ad,Ao,icol,ierr))
 31:       PetscCallA(MatGetOwnershipRange(A,rstart,rend,ierr))
 32: !
 33: !   Print diagonal portion of matrix
 34: !
 35:       PetscCallA(PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1,ierr))
 36:       zero = 0
 37:       one  = 1
 38:       bb = .true.
 39:       PetscCallA(MatGetRowIJF90(Ad,one,bb,bb,n,ia,ja,done,ierr))
 40:       PetscCallA(MatSeqAIJGetArrayF90(Ad,aa,ierr))
 41:       do 10, i=1,n
 42:         write(7+rank,*) 'row ',i+rstart,' number nonzeros ',ia(i+1)-ia(i)
 43:         do 20, j=ia(i),ia(i+1)-1
 44:           write(7+rank,*)'  ',j,ja(j)+rstart,aa(j)
 45:  20     continue
 46:  10   continue
 47:       PetscCallA(MatRestoreRowIJF90(Ad,one,bb,bb,n,ia,ja,done,ierr))
 48:       PetscCallA(MatSeqAIJRestoreArrayF90(Ad,aa,ierr))
 49: !
 50: !   Print off-diagonal portion of matrix
 51: !
 52:       PetscCallA(MatGetRowIJF90(Ao,one,bb,bb,n,ia,ja,done,ierr))
 53:       PetscCallA(MatSeqAIJGetArrayF90(Ao,aa,ierr))
 54:       do 30, i=1,n
 55:         write(7+rank,*) 'row ',i+rstart,' number nonzeros ',ia(i+1)-ia(i)
 56:         do 40, j=ia(i),ia(i+1)-1
 57:           write(7+rank,*)'  ',j,icol(ja(j))+1,aa(j)
 58:  40     continue
 59:  30   continue
 60:       PetscCallA(MatMPIAIJRestoreSeqAIJF90(A,Ad,Ao,icol,ierr))
 61:       PetscCallA(MatRestoreRowIJF90(Ao,one,bb,bb,n,ia,ja,done,ierr))
 62:       PetscCallA(MatSeqAIJRestoreArrayF90(Ao,aa,ierr))

 64:       PetscCallA(PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1,ierr))

 66:       PetscCallA(MatGetDiagonalBlock(A,Ad,ierr))
 67:       PetscCallA(MatView(Ad,PETSC_VIEWER_STDOUT_WORLD,ierr))

 69:       PetscCallA(MatView(A,PETSC_VIEWER_STDOUT_WORLD,ierr))
 70:       PetscCallA(MatDestroy(A,ierr))
 71:       PetscCallA(PetscViewerDestroy(v,ierr))
 72:       PetscCallA(PetscFinalize(ierr))
 73:       end

 75: !/*TEST
 76: !
 77: !     build:
 78: !       requires: double !complex !defined(PETSC_USE_64BIT_INDICES)
 79: !
 80: !     test:
 81: !        args: -binary_read_double -options_left false
 82: !
 83: !TEST*/