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