1*d8606c27SBarry Smith! 2*d8606c27SBarry Smith! This program demonstrates use of MatCreateSubMatrices() from Fortran 3*d8606c27SBarry Smith! 4*d8606c27SBarry Smith program main 5*d8606c27SBarry Smith#include <petsc/finclude/petscmat.h> 6*d8606c27SBarry Smith use petscmat 7*d8606c27SBarry Smith implicit none 8*d8606c27SBarry Smith 9*d8606c27SBarry Smith Mat A,B(2) 10*d8606c27SBarry Smith PetscErrorCode ierr 11*d8606c27SBarry Smith PetscInt nis,zero(1) 12*d8606c27SBarry Smith PetscViewer v 13*d8606c27SBarry Smith IS isrow 14*d8606c27SBarry Smith PetscMPIInt rank 15*d8606c27SBarry Smith 16*d8606c27SBarry Smith PetscCallA(PetscInitialize(ierr)) 17*d8606c27SBarry Smith PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)) 18*d8606c27SBarry Smith 19*d8606c27SBarry Smith#if defined(PETSC_USE_64BIT_INDICES) 20*d8606c27SBarry Smith PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD,'${PETSC_DIR}/share/petsc/datafiles/matrices/' // 'ns-real-int64-float64',FILE_MODE_READ,v,ierr)) 21*d8606c27SBarry Smith#else 22*d8606c27SBarry Smith PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD,'${PETSC_DIR}/share/petsc/datafiles/matrices/' // 'ns-real-int32-float64',FILE_MODE_READ,v,ierr)) 23*d8606c27SBarry Smith#endif 24*d8606c27SBarry Smith 25*d8606c27SBarry Smith PetscCallA(MatCreate(PETSC_COMM_WORLD,A,ierr)) 26*d8606c27SBarry Smith PetscCallA(MatSetType(A, MATAIJ,ierr)) 27*d8606c27SBarry Smith PetscCallA(MatLoad(A,v,ierr)) 28*d8606c27SBarry Smith 29*d8606c27SBarry Smith nis = 1 30*d8606c27SBarry Smith zero(1) = 0 31*d8606c27SBarry Smith if (rank .eq. 1) then 32*d8606c27SBarry Smith nis = 0 ! test nis = 0 33*d8606c27SBarry Smith endif 34*d8606c27SBarry Smith PetscCallA(ISCreateGeneral(PETSC_COMM_SELF,nis,zero,PETSC_COPY_VALUES,isrow,ierr)) 35*d8606c27SBarry Smith 36*d8606c27SBarry Smith PetscCallA(MatCreateSubmatrices(A,nis,isrow,isrow,MAT_INITIAL_MATRIX,B,ierr)) 37*d8606c27SBarry Smith 38*d8606c27SBarry Smith if (rank .eq. 0) then 39*d8606c27SBarry Smith PetscCallA(MatView(B(1),PETSC_VIEWER_STDOUT_SELF,ierr)) 40*d8606c27SBarry Smith endif 41*d8606c27SBarry Smith 42*d8606c27SBarry Smith PetscCallA(MatCreateSubmatrices(A,nis,isrow,isrow,MAT_REUSE_MATRIX,B,ierr)) 43*d8606c27SBarry Smith 44*d8606c27SBarry Smith if (rank .eq. 0) then 45*d8606c27SBarry Smith PetscCallA(MatView(B(1),PETSC_VIEWER_STDOUT_SELF,ierr)) 46*d8606c27SBarry Smith endif 47*d8606c27SBarry Smith 48*d8606c27SBarry Smith PetscCallA(ISDestroy(isrow,ierr)) 49*d8606c27SBarry Smith PetscCallA(MatDestroy(A,ierr)) 50*d8606c27SBarry Smith PetscCallA(MatDestroySubMatrices(nis,B,ierr)) 51*d8606c27SBarry Smith PetscCallA(PetscViewerDestroy(v,ierr)) 52*d8606c27SBarry Smith 53*d8606c27SBarry Smith PetscCallA(PetscFinalize(ierr)) 54*d8606c27SBarry Smith end 55*d8606c27SBarry Smith 56*d8606c27SBarry Smith!/*TEST 57*d8606c27SBarry Smith! 58*d8606c27SBarry Smith! test: 59*d8606c27SBarry Smith! requires: double !complex 60*d8606c27SBarry Smith! 61*d8606c27SBarry Smith!TEST*/ 62*d8606c27SBarry Smith 63