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