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