1! 2! This program demonstrates use of MatGetRowIJ() from Fortran 3! 4#include <petsc/finclude/petscmat.h> 5program main 6 7 use petscmat 8 implicit none 9 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(:) 20 21 PetscCallA(PetscInitialize(ierr)) 22 PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr)) 23 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)) 29 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 i = 1, n 42 write (7 + rank, *) 'row ', i + rstart, ' number nonzeros ', ia(i + 1) - ia(i) 43 do j = ia(i), ia(i + 1) - 1 44 write (7 + rank, *) ' ', j, ja(j) + rstart, aa(j) 45 end do 46 end do 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 i = 1, n 55 write (7 + rank, *) 'row ', i + rstart, ' number nonzeros ', ia(i + 1) - ia(i) 56 do j = ia(i), ia(i + 1) - 1 57 write (7 + rank, *) ' ', j, icol(ja(j)) + 1, aa(j) 58 end do 59 end do 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)) 63 64 PetscCallA(PetscSequentialPhaseEnd(PETSC_COMM_WORLD, 1, ierr)) 65 66 PetscCallA(MatGetDiagonalBlock(A, Ad, ierr)) 67 PetscCallA(MatView(Ad, PETSC_VIEWER_STDOUT_WORLD, ierr)) 68 69 PetscCallA(MatView(A, PETSC_VIEWER_STDOUT_WORLD, ierr)) 70 PetscCallA(MatDestroy(A, ierr)) 71 PetscCallA(PetscViewerDestroy(v, ierr)) 72 PetscCallA(PetscFinalize(ierr)) 73end 74 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*/ 84