Actual source code: ex37.c
petsc-3.14.6 2021-03-30
2: static char help[] = "Tests MatCopy() and MatStore/RetrieveValues().\n\n";
4: #include <petscmat.h>
6: int main(int argc,char **args)
7: {
8: Mat C,A;
9: PetscInt i, n = 10,midx[3],bs=1;
11: PetscScalar v[3];
12: PetscBool flg,isAIJ;
13: MatType type;
14: PetscMPIInt size;
16: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
17: MPI_Comm_size(PETSC_COMM_WORLD,&size);
18: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
19: PetscOptionsGetInt(NULL,NULL,"-mat_block_size",&bs,NULL);
21: MatCreate(PETSC_COMM_WORLD,&C);
22: MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,n,n);
23: MatSetType(C,MATAIJ);
24: MatSetFromOptions(C);
25: PetscObjectSetName((PetscObject)C,"initial");
27: MatGetType(C,&type);
28: if (size == 1) {
29: PetscObjectTypeCompare((PetscObject)C,MATSEQAIJ,&isAIJ);
30: } else {
31: PetscObjectTypeCompare((PetscObject)C,MATMPIAIJ,&isAIJ);
32: }
33: MatSeqAIJSetPreallocation(C,3,NULL);
34: MatMPIAIJSetPreallocation(C,3,NULL,3,NULL);
35: MatSeqBAIJSetPreallocation(C,bs,3,NULL);
36: MatMPIBAIJSetPreallocation(C,bs,3,NULL,3,NULL);
37: MatSeqSBAIJSetPreallocation(C,bs,3,NULL);
38: MatMPISBAIJSetPreallocation(C,bs,3,NULL,3,NULL);
40: v[0] = -1.; v[1] = 2.; v[2] = -1.;
41: for (i=1; i<n-1; i++) {
42: midx[2] = i-1; midx[1] = i; midx[0] = i+1;
43: MatSetValues(C,1,&i,3,midx,v,INSERT_VALUES);
44: }
45: i = 0; midx[0] = 0; midx[1] = 1;
46: v[0] = 2.0; v[1] = -1.;
47: MatSetValues(C,1,&i,2,midx,v,INSERT_VALUES);
48: i = n-1; midx[0] = n-2; midx[1] = n-1;
49: v[0] = -1.0; v[1] = 2.;
50: MatSetValues(C,1,&i,2,midx,v,INSERT_VALUES);
52: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
53: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
54: MatView(C,NULL);
55: MatViewFromOptions(C,NULL,"-view");
57: /* test matduplicate */
58: MatDuplicate(C,MAT_COPY_VALUES,&A);
59: PetscObjectSetName((PetscObject)A,"duplicate_copy");
60: MatViewFromOptions(A,NULL,"-view");
61: MatEqual(A,C,&flg);
62: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"MatDuplicate(C,MAT_COPY_VALUES,): Matrices are NOT equal");
63: MatDestroy(&A);
65: /* test matrices with different nonzero patterns - Note: A is created with different nonzero pattern of C! */
66: MatCreate(PETSC_COMM_WORLD,&A);
67: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
68: MatSetFromOptions(A);
69: MatSetUp(A);
70: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
71: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
73: MatCopy(C,A,DIFFERENT_NONZERO_PATTERN);
74: PetscObjectSetName((PetscObject)A,"copy_diffnnz");
75: MatViewFromOptions(A,NULL,"-view");
76: MatEqual(A,C,&flg);
77: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"MatCopy(C,A,DIFFERENT_NONZERO_PATTERN): Matrices are NOT equal");
79: /* test matrices with same nonzero pattern */
80: MatDestroy(&A);
81: MatDuplicate(C,MAT_DO_NOT_COPY_VALUES,&A);
82: MatCopy(C,A,SAME_NONZERO_PATTERN);
83: PetscObjectSetName((PetscObject)A,"copy_samennz");
84: MatViewFromOptions(A,NULL,"-view");
85: MatEqual(A,C,&flg);
86: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"MatCopy(C,A,SAME_NONZERO_PATTERN): Matrices are NOT equal");
88: /* test subset nonzero pattern */
89: MatCopy(C,A,SUBSET_NONZERO_PATTERN);
90: PetscObjectSetName((PetscObject)A,"copy_subnnz");
91: MatViewFromOptions(A,NULL,"-view");
92: MatEqual(A,C,&flg);
93: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"MatCopy(C,A,SUBSET_NONZERO_PATTERN): Matrices are NOT equal");
95: /* Test MatCopy on a matrix obtained after MatConvert from AIJ
96: see https://lists.mcs.anl.gov/pipermail/petsc-dev/2019-April/024289.html */
97: MatHasCongruentLayouts(C,&flg);
98: if (flg) {
99: Mat Cs,Cse;
100: MatType Ctype,Cstype;
102: MatGetType(C,&Ctype);
103: MatTranspose(C,MAT_INITIAL_MATRIX,&Cs);
104: MatAXPY(Cs,1.0,C,DIFFERENT_NONZERO_PATTERN);
105: MatConvert(Cs,MATAIJ,MAT_INPLACE_MATRIX,&Cs);
106: MatSetOption(Cs,MAT_SYMMETRIC,PETSC_TRUE);
107: MatGetType(Cs,&Cstype);
109: PetscObjectSetName((PetscObject)Cs,"symm_initial");
110: MatViewFromOptions(Cs,NULL,"-view");
112: MatConvert(Cs,Ctype,MAT_INITIAL_MATRIX,&Cse);
113: PetscObjectSetName((PetscObject)Cse,"symm_conv_init");
114: MatViewFromOptions(Cse,NULL,"-view");
115: MatMultEqual(Cs,Cse,5,&flg);
116: if (!flg) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatConvert MAT_INITIAL_MATRIX %s -> %s: Matrices are NOT multequal",Ctype,Cstype);
118: MatConvert(Cs,Ctype,MAT_REUSE_MATRIX,&Cse);
119: PetscObjectSetName((PetscObject)Cse,"symm_conv_reuse");
120: MatViewFromOptions(Cse,NULL,"-view");
121: MatMultEqual(Cs,Cse,5,&flg);
122: if (!flg) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatConvert MAT_REUSE_MATRIX %s -> %s: Matrices are NOT multequal",Ctype,Cstype);
124: MatCopy(Cs,Cse,SAME_NONZERO_PATTERN);
125: PetscObjectSetName((PetscObject)Cse,"symm_conv_copy_samennz");
126: MatViewFromOptions(Cse,NULL,"-view");
127: MatMultEqual(Cs,Cse,5,&flg);
128: if (!flg) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatCopy(...SAME_NONZERO_PATTERN) %s -> %s: Matrices are NOT multequal",Ctype,Cstype);
130: MatCopy(Cs,Cse,SUBSET_NONZERO_PATTERN);
131: PetscObjectSetName((PetscObject)Cse,"symm_conv_copy_subnnz");
132: MatViewFromOptions(Cse,NULL,"-view");
133: MatMultEqual(Cs,Cse,5,&flg);
134: if (!flg) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatCopy(...SUBSET_NONZERO_PATTERN) %s -> %s: Matrices are NOT multequal",Ctype,Cstype);
136: MatCopy(Cs,Cse,DIFFERENT_NONZERO_PATTERN);
137: PetscObjectSetName((PetscObject)Cse,"symm_conv_copy_diffnnz");
138: MatViewFromOptions(Cse,NULL,"-view");
139: MatMultEqual(Cs,Cse,5,&flg);
140: if (!flg) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatCopy(...DIFFERENT_NONZERO_PATTERN) %s -> %s: Matrices are NOT multequal",Ctype,Cstype);
142: MatDestroy(&Cse);
143: MatDestroy(&Cs);
144: }
146: /* test MatStore/RetrieveValues() */
147: if (isAIJ) {
148: MatSetOption(A,MAT_NEW_NONZERO_LOCATIONS,PETSC_FALSE);
149: MatStoreValues(A);
150: MatZeroEntries(A);
151: MatRetrieveValues(A);
152: }
154: MatDestroy(&C);
155: MatDestroy(&A);
156: PetscFinalize();
157: return ierr;
158: }
163: /*TEST
165: testset:
166: nsize: {{1 2}separate output}
167: args: -view ::ascii_info -mat_type {{aij baij sbaij mpiaij mpibaij mpisbaij}separate output} -mat_block_size {{1 2}separate output}
169: TEST*/