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