Actual source code: ex90.c
2: static char help[] = "Tests MatPtAP() \n";
4: #include <petscmat.h>
6: /*
7: * This is an extremely simple example to test MatPtAP. It is very useful when developing and debugging the code.
8: *
9: * A =
11: 1 2 0 4
12: 0 1 2 0
13: 2 0 4 0
14: 0 1 2 1
15: *
16: *
17: *
18: * P =
20: 1.00000 0.00000
21: 0.30000 0.50000
22: 0.00000 0.80000
23: 0.90000 0.00000
24: *
25: *
26: *AP =
27: *
28: * 5.20000 1.00000
29: 0.30000 2.10000
30: 2.00000 3.20000
31: 1.20000 2.10000
32: *
33: * PT =
35: 1.00000 0.30000 0.00000 0.90000
36: 0.00000 0.50000 0.80000 0.00000
38: *
39: * C =
41: 6.3700 3.5200
42: 1.7500 3.6100
43: *
44: * */
46: int main(int argc, char **argv)
47: {
48: Mat A, P, PtAP;
49: PetscInt i1[] = {0, 3, 5}, i2[] = {0, 2, 5};
50: PetscInt j1[] = {0, 1, 3, 1, 2}, j2[] = {0, 2, 1, 2, 3};
51: PetscScalar a1[] = {1, 2, 4, 1, 2}, a2[] = {2, 4, 1, 2, 1};
52: PetscInt pi1[] = {0, 1, 3}, pi2[] = {0, 1, 2};
53: PetscInt pj1[] = {0, 0, 1}, pj2[] = {1, 0};
54: PetscScalar pa1[] = {1, 0.3, 0.5}, pa2[] = {0.8, 0.9};
55: MPI_Comm comm;
56: PetscMPIInt rank, size;
59: PetscInitialize(&argc, &argv, NULL, help);
60: comm = PETSC_COMM_WORLD;
61: MPI_Comm_rank(comm, &rank);
62: MPI_Comm_size(comm, &size);
64: MatCreateMPIAIJWithArrays(comm, 2, 2, PETSC_DETERMINE, PETSC_DETERMINE, rank ? i2 : i1, rank ? j2 : j1, rank ? a2 : a1, &A);
65: MatCreateMPIAIJWithArrays(comm, 2, 1, PETSC_DETERMINE, PETSC_DETERMINE, rank ? pi2 : pi1, rank ? pj2 : pj1, rank ? pa2 : pa1, &P);
66: MatPtAP(A, P, MAT_INITIAL_MATRIX, 1.1, &PtAP);
67: MatView(A, NULL);
68: MatView(P, NULL);
69: MatView(PtAP, NULL);
70: MatPtAP(A, P, MAT_REUSE_MATRIX, 1.1, &PtAP);
71: MatView(A, NULL);
72: MatView(P, NULL);
73: MatView(PtAP, NULL);
74: MatDestroy(&A);
75: MatDestroy(&P);
76: MatDestroy(&PtAP);
77: PetscFinalize();
78: return 0;
79: }
81: /*TEST
82: test:
83: nsize: 2
84: args: -matptap_via allatonce
85: output_file: output/ex90_1.out
87: test:
88: nsize: 2
89: suffix: merged
90: args: -matptap_via allatonce_merged
91: output_file: output/ex90_1.out
93: TEST*/