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