Actual source code: ex79f.F90
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(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))
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*/