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(MatMPIAIJGetSeqAIJ(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(MatGetRowIJ(Ad, one, bb, bb, n, ia, ja, done, ierr))
 40:   PetscCallA(MatSeqAIJGetArray(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(MatRestoreRowIJ(Ad, one, bb, bb, n, ia, ja, done, ierr))
 48:       PetscCallA(MatSeqAIJRestoreArray(Ad, aa, ierr))
 49: !
 50: !   Print off-diagonal portion of matrix
 51: !
 52:       PetscCallA(MatGetRowIJ(Ao, one, bb, bb, n, ia, ja, done, ierr))
 53:       PetscCallA(MatSeqAIJGetArray(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(MatMPIAIJRestoreSeqAIJ(A, Ad, Ao, icol, ierr))
 61:           PetscCallA(MatRestoreRowIJ(Ao, one, bb, bb, n, ia, ja, done, ierr))
 62:           PetscCallA(MatSeqAIJRestoreArray(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*/