Actual source code: ex95.c
1: static char help[] = "Testing MatCreateMPIAIJSumSeqAIJ().\n\n";
3: #include <petscmat.h>
5: int main(int argc, char **argv)
6: {
7: Mat A, B;
8: MatScalar a[1], alpha;
9: PetscMPIInt size, rank;
10: PetscInt m, n, i, col, prid;
13: PetscInitialize(&argc, &argv, (char *)0, help);
14: MPI_Comm_size(PETSC_COMM_WORLD, &size);
15: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
16: prid = size;
17: PetscOptionsGetInt(NULL, NULL, "-prid", &prid, NULL);
19: m = n = 10 * size;
20: MatCreate(PETSC_COMM_SELF, &A);
21: MatSetSizes(A, PETSC_DETERMINE, PETSC_DETERMINE, m, n);
22: MatSetType(A, MATSEQAIJ);
23: MatSetUp(A);
25: a[0] = rank + 1;
26: for (i = 0; i < m - rank; i++) {
27: col = i + rank;
28: MatSetValues(A, 1, &i, 1, &col, a, INSERT_VALUES);
29: }
30: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
31: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
33: if (rank == prid) {
34: PetscPrintf(PETSC_COMM_SELF, "[%d] A: \n", rank);
35: MatView(A, PETSC_VIEWER_STDOUT_SELF);
36: }
38: /* Test MatCreateMPIAIJSumSeqAIJ */
39: MatCreateMPIAIJSumSeqAIJ(PETSC_COMM_WORLD, A, PETSC_DECIDE, PETSC_DECIDE, MAT_INITIAL_MATRIX, &B);
41: /* Test MAT_REUSE_MATRIX */
42: alpha = 0.1;
43: for (i = 0; i < 3; i++) {
44: MatScale(A, alpha);
45: MatCreateMPIAIJSumSeqAIJ(PETSC_COMM_WORLD, A, PETSC_DECIDE, PETSC_DECIDE, MAT_REUSE_MATRIX, &B);
46: }
47: MatView(B, PETSC_VIEWER_STDOUT_WORLD);
48: MatDestroy(&B);
49: MatDestroy(&A);
50: PetscFinalize();
51: return 0;
52: }
54: /*TEST
56: test:
57: nsize: 3
58: filter: grep -v " MPI process"
60: test:
61: suffix: 2
62: filter: grep -v " MPI process"
64: TEST*/