1! 2! This program demonstrates use of MatGetRowIJ() from Fortran 3! 4#include <petsc/finclude/petscmat.h> 5program main 6 use petscmat 7 implicit none 8 9 Mat A, Ad, Ao 10 PetscErrorCode ierr 11 PetscMPIInt rank 12 PetscViewer v 13 PetscInt i, j 14 PetscInt n, rstart 15 PetscInt zero, one, rend 16 PetscBool done, bb 17 PetscScalar, pointer :: aa(:) 18 PetscInt, pointer :: ia(:), ja(:), icol(:) 19 20 PetscCallA(PetscInitialize(ierr)) 21 PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr)) 22 23 PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD, '${PETSC_DIR}/share/petsc/datafiles/matrices/'//'ns-real-int32-float64', FILE_MODE_READ, v, ierr)) 24 PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr)) 25 PetscCallA(MatSetType(A, MATMPIAIJ, ierr)) 26 PetscCallA(MatLoad(A, v, ierr)) 27 PetscCallA(MatView(A, PETSC_VIEWER_STDOUT_WORLD, ierr)) 28 29 PetscCallA(MatMPIAIJGetSeqAIJ(A, Ad, Ao, icol, ierr)) 30 PetscCallA(MatGetOwnershipRange(A, rstart, rend, ierr)) 31! 32! Print diagonal portion of matrix 33! 34 PetscCallA(PetscSequentialPhaseBegin(PETSC_COMM_WORLD, 1, ierr)) 35 zero = 0 36 one = 1 37 bb = .true. 38 PetscCallA(MatGetRowIJ(Ad, one, bb, bb, n, ia, ja, done, ierr)) 39 PetscCallA(MatSeqAIJGetArray(Ad, aa, ierr)) 40 do i = 1, n 41 write (7 + rank, *) 'row ', i + rstart, ' number nonzeros ', ia(i + 1) - ia(i) 42 do j = ia(i), ia(i + 1) - 1 43 write (7 + rank, *) ' ', j, ja(j) + rstart, aa(j) 44 end do 45 end do 46 PetscCallA(MatRestoreRowIJ(Ad, one, bb, bb, n, ia, ja, done, ierr)) 47 PetscCallA(MatSeqAIJRestoreArray(Ad, aa, ierr)) 48! 49! Print off-diagonal portion of matrix 50! 51 PetscCallA(MatGetRowIJ(Ao, one, bb, bb, n, ia, ja, done, ierr)) 52 PetscCallA(MatSeqAIJGetArray(Ao, aa, ierr)) 53 do i = 1, n 54 write (7 + rank, *) 'row ', i + rstart, ' number nonzeros ', ia(i + 1) - ia(i) 55 do j = ia(i), ia(i + 1) - 1 56 write (7 + rank, *) ' ', j, icol(ja(j)) + 1, aa(j) 57 end do 58 end do 59 PetscCallA(MatMPIAIJRestoreSeqAIJ(A, Ad, Ao, icol, ierr)) 60 PetscCallA(MatRestoreRowIJ(Ao, one, bb, bb, n, ia, ja, done, ierr)) 61 PetscCallA(MatSeqAIJRestoreArray(Ao, aa, ierr)) 62 63 PetscCallA(PetscSequentialPhaseEnd(PETSC_COMM_WORLD, 1, ierr)) 64 65 PetscCallA(MatGetDiagonalBlock(A, Ad, ierr)) 66 PetscCallA(MatView(Ad, PETSC_VIEWER_STDOUT_WORLD, ierr)) 67 68 PetscCallA(MatView(A, PETSC_VIEWER_STDOUT_WORLD, ierr)) 69 PetscCallA(MatDestroy(A, ierr)) 70 PetscCallA(PetscViewerDestroy(v, ierr)) 71 PetscCallA(PetscFinalize(ierr)) 72end 73 74!/*TEST 75! 76! build: 77! requires: double !complex !defined(PETSC_USE_64BIT_INDICES) 78! 79! test: 80! args: -binary_read_double -options_left false 81! 82!TEST*/ 83