Actual source code: ex102.c
2: static char help[] = "Tests MatCreateLRC()\n\n";
4: #include <petscmat.h>
6: int main(int argc, char **args)
7: {
8: Vec x, b, c = NULL;
9: Mat A, U, V, LR, X, LRe;
10: PetscInt M = 5, N = 7;
11: PetscBool flg;
14: PetscInitialize(&argc, &args, (char *)0, help);
15: PetscOptionsGetInt(NULL, NULL, "-m", &M, NULL);
16: PetscOptionsGetInt(NULL, NULL, "-n", &N, NULL);
18: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
19: Create the sparse matrix
20: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
21: MatCreate(PETSC_COMM_WORLD, &A);
22: MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, M, N);
23: MatSetOptionsPrefix(A, "A_");
24: MatSetFromOptions(A);
25: MatSeqAIJSetPreallocation(A, 5, NULL);
26: MatMPIAIJSetPreallocation(A, 5, NULL, 5, NULL);
27: MatSetUp(A);
28: MatSetRandom(A, NULL);
30: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
31: Create the dense matrices
32: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
33: MatCreate(PETSC_COMM_WORLD, &U);
34: MatSetSizes(U, PETSC_DECIDE, PETSC_DECIDE, M, 3);
35: MatSetType(U, MATDENSE);
36: MatSetOptionsPrefix(U, "U_");
37: MatSetFromOptions(U);
38: MatSetUp(U);
39: MatSetRandom(U, NULL);
41: MatCreate(PETSC_COMM_WORLD, &V);
42: MatSetSizes(V, PETSC_DECIDE, PETSC_DECIDE, N, 3);
43: MatSetType(V, MATDENSE);
44: MatSetOptionsPrefix(V, "V_");
45: MatSetFromOptions(V);
46: MatSetUp(V);
47: MatSetRandom(V, NULL);
49: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
50: Create a vector to hold the diagonal of C
51: A sequential vector can be created as well on each process
52: It is user responsibility to ensure the data in the vector
53: is consistent across processors
54: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
55: PetscOptionsHasName(NULL, NULL, "-use_c", &flg);
56: if (flg) {
57: VecCreateMPI(PETSC_COMM_WORLD, PETSC_DECIDE, 3, &c);
58: VecSetRandom(c, NULL);
59: }
61: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
62: Create low rank correction matrix
63: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
64: PetscOptionsHasName(NULL, NULL, "-low_rank", &flg);
65: if (flg) {
66: /* create a low-rank matrix, with no A-matrix */
67: MatCreateLRC(NULL, U, c, V, &LR);
68: MatDestroy(&A);
69: } else {
70: MatCreateLRC(A, U, c, V, &LR);
71: }
73: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: Create the low rank correction matrix explicitly to check for
75: correctness
76: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
77: MatHermitianTranspose(V, MAT_INITIAL_MATRIX, &X);
78: MatDiagonalScale(X, c, NULL);
79: MatMatMult(U, X, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &LRe);
80: MatDestroy(&X);
81: if (A) MatAYPX(LRe, 1.0, A, DIFFERENT_NONZERO_PATTERN);
83: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
84: Create test vectors
85: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
86: MatCreateVecs(LR, &x, &b);
87: VecSetRandom(x, NULL);
88: MatMult(LR, x, b);
89: MatMultTranspose(LR, b, x);
90: VecDestroy(&x);
91: VecDestroy(&b);
93: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
94: Check correctness
95: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
96: MatMultEqual(LR, LRe, 10, &flg);
97: if (!flg) PetscPrintf(PETSC_COMM_WORLD, "Error in MatMult\n");
98: #if !defined(PETSC_USE_COMPLEX)
99: MatMultHermitianTransposeEqual(LR, LRe, 10, &flg);
100: if (!flg) PetscPrintf(PETSC_COMM_WORLD, "Error in MatMultTranspose\n");
101: #endif
103: MatDestroy(&A);
104: MatDestroy(&LRe);
105: MatDestroy(&U);
106: MatDestroy(&V);
107: VecDestroy(&c);
108: MatDestroy(&LR);
110: /*
111: Always call PetscFinalize() before exiting a program. This routine
112: - finalizes the PETSc libraries as well as MPI
113: - provides summary and diagnostic information if certain runtime
114: options are chosen (e.g., -log_view).
115: */
116: PetscFinalize();
117: return 0;
118: }
120: /*TEST
122: testset:
123: output_file: output/ex102_1.out
124: nsize: {{1 2}}
125: args: -low_rank {{0 1}} -use_c {{0 1}}
126: test:
127: suffix: standard
128: test:
129: suffix: cuda
130: requires: cuda
131: args: -A_mat_type aijcusparse -U_mat_type densecuda -V_mat_type densecuda
133: TEST*/