1! 2! This program demonstrates use of MatGetRowIJ() from Fortran 3! 4 program main 5 6#include <petsc/finclude/petscmat.h> 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(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)) 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)) 73 end 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