xref: /petsc/src/mat/tests/ex249.c (revision ebead697dbf761eb322f829370bbe90b3bd93fa3)
1 static char help[] = "Test MatCreateSubMatrices\n\n";
2 
3 #include <petscis.h>
4 #include <petscmat.h>
5 
6 int main(int argc,char **args)
7 {
8   Mat             A,*submats,*submats2;
9   IS              *irow,*icol;
10   PetscInt        i,n;
11   PetscMPIInt     rank;
12   PetscViewer     matfd,rowfd,colfd;
13   PetscBool       same;
14   char            matfile[PETSC_MAX_PATH_LEN],rowfile[PETSC_MAX_PATH_LEN],colfile[PETSC_MAX_PATH_LEN];
15   char            rankstr[16]={0};
16 
17   PetscFunctionBeginUser;
18   PetscCall(PetscInitialize(&argc,&args,(char*)0,help));
19   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank));
20 
21   PetscCall(PetscOptionsGetString(NULL,NULL,"-A",matfile,sizeof(matfile),NULL));
22   PetscCall(PetscOptionsGetString(NULL,NULL,"-row",rowfile,sizeof(rowfile),NULL));
23   PetscCall(PetscOptionsGetString(NULL,NULL,"-col",colfile,sizeof(colfile),NULL));
24 
25   /* Each rank has its own files for row/col ISes */
26   PetscCall(PetscSNPrintf(rankstr,16,"-%d",rank));
27   PetscCall(PetscStrlcat(rowfile,rankstr,PETSC_MAX_PATH_LEN));
28   PetscCall(PetscStrlcat(colfile,rankstr,PETSC_MAX_PATH_LEN));
29 
30   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,matfile,FILE_MODE_READ,&matfd));
31   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF,rowfile,FILE_MODE_READ,&rowfd));
32   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF,colfile,FILE_MODE_READ,&colfd));
33 
34   PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
35   PetscCall(MatSetFromOptions(A));
36   PetscCall(MatLoad(A,matfd));
37 
38   /* We stored the number of ISes at the beginning of rowfd */
39   PetscCall(PetscViewerBinaryRead(rowfd,&n,1,NULL,PETSC_INT));
40   PetscCall(PetscMalloc2(n,&irow,n,&icol));
41   for (i=0; i<n; i++) {
42     PetscCall(ISCreate(PETSC_COMM_SELF,&irow[i]));
43     PetscCall(ISCreate(PETSC_COMM_SELF,&icol[i]));
44     PetscCall(ISLoad(irow[i],rowfd));
45     PetscCall(ISLoad(icol[i],colfd));
46   }
47 
48   PetscCall(PetscViewerDestroy(&matfd));
49   PetscCall(PetscViewerDestroy(&rowfd));
50   PetscCall(PetscViewerDestroy(&colfd));
51 
52   /* Create submats for the first time */
53   PetscCall(MatCreateSubMatrices(A,n,irow,icol,MAT_INITIAL_MATRIX,&submats));
54 
55   /* Dup submats to submats2 for later comparison */
56   PetscCall(PetscMalloc1(n,&submats2));
57   for (i=0; i<n; i++) {
58     PetscCall(MatDuplicate(submats[i],MAT_COPY_VALUES,&submats2[i]));
59   }
60 
61   /* Create submats again */
62   PetscCall(MatCreateSubMatrices(A,n,irow,icol,MAT_REUSE_MATRIX,&submats));
63 
64   /* Compare submats and submats2 */
65   for (i=0; i<n; i++) {
66     PetscCall(MatEqual(submats[i],submats2[i],&same));
67     PetscCheck(same,PETSC_COMM_SELF,PETSC_ERR_PLIB,"submatrix %" PetscInt_FMT " is not same",i);
68   }
69 
70   PetscCall(MatDestroy(&A));
71   for (i=0; i<n; i++) {
72     PetscCall(ISDestroy(&irow[i]));
73     PetscCall(ISDestroy(&icol[i]));
74   }
75   PetscCall(MatDestroySubMatrices(n,&submats));
76   PetscCall(MatDestroyMatrices(n,&submats2));
77   PetscCall(PetscFree2(irow,icol));
78   PetscCall(PetscFinalize());
79   return 0;
80 }
81 
82 /*TEST
83 
84    test:
85      suffix: 1
86      nsize: 2
87      requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
88      args: -mat_type {{aij baij}} -A ${DATAFILESPATH}/matrices/CreateSubMatrices/A -row ${DATAFILESPATH}/matrices/CreateSubMatrices/row -col ${DATAFILESPATH}/matrices/CreateSubMatrices/col
89 
90 TEST*/
91