Actual source code: ex197.c
1: static char help[] = "Test MatMultHermitianTranspose() and MatMultHermitianTransposeAdd().\n\n";
3: #include <petscmat.h>
5: int main(int argc, char **args)
6: {
7: Mat A, B, C;
8: Vec x, y, ys;
9: PetscInt i, j;
10: PetscScalar v;
11: PetscBool flg;
13: PetscFunctionBeginUser;
14: PetscCall(PetscInitialize(&argc, &args, NULL, help));
15: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
16: PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, 2, 2));
17: PetscCall(MatSetType(A, MATAIJ));
18: PetscCall(MatSetFromOptions(A));
19: PetscCall(MatSetUp(A));
21: i = 0;
22: j = 0;
23: v = 2.0;
24: PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES));
25: i = 0;
26: j = 1;
27: v = 3.0 + 4.0 * PETSC_i;
28: PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES));
29: i = 1;
30: j = 0;
31: v = 5.0 + 6.0 * PETSC_i;
32: PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES));
33: i = 1;
34: j = 1;
35: v = 7.0 + 8.0 * PETSC_i;
36: PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES));
37: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
38: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
40: /* Create vectors */
41: PetscCall(MatCreateVecs(A, &x, &y));
42: PetscCall(VecDuplicate(y, &ys));
44: i = 0;
45: v = 10.0 + 11.0 * PETSC_i;
46: PetscCall(VecSetValues(x, 1, &i, &v, INSERT_VALUES));
47: i = 1;
48: v = 100.0 + 120.0 * PETSC_i;
49: PetscCall(VecSetValues(x, 1, &i, &v, INSERT_VALUES));
50: PetscCall(VecAssemblyBegin(x));
51: PetscCall(VecAssemblyEnd(x));
53: PetscCall(MatMultHermitianTranspose(A, x, y));
54: PetscCall(VecView(y, PETSC_VIEWER_STDOUT_WORLD));
55: PetscCall(MatMultHermitianTransposeAdd(A, x, y, ys));
56: PetscCall(VecView(ys, PETSC_VIEWER_STDOUT_WORLD));
58: PetscCall(MatHermitianTranspose(A, MAT_INITIAL_MATRIX, &B));
59: PetscCall(MatCreateHermitianTranspose(A, &C));
60: PetscCall(MatMultHermitianTransposeEqual(B, C, 4, &flg));
61: PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "B^Hx != C^Hx");
62: PetscCall(MatMultHermitianTransposeAddEqual(B, C, 4, &flg));
63: PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "y+B^Hx != y+C^Hx");
64: PetscCall(MatDestroy(&C));
65: PetscCall(MatDestroy(&B));
67: PetscCall(MatDestroy(&A));
69: PetscCall(VecDestroy(&x));
70: PetscCall(VecDestroy(&y));
71: PetscCall(VecDestroy(&ys));
72: PetscCall(PetscFinalize());
73: return 0;
74: }
76: /*TEST
78: build:
79: requires: complex
81: testset:
82: output_file: output/ex197_1.out
83: test:
84: suffix: 1
85: args: -mat_type {{aij dense}}
86: test:
87: suffix: 1_cuda
88: requires: cuda
89: args: -mat_type densecuda
90: filter: sed -e 's/seqcuda/seq/'
92: testset:
93: output_file: output/ex197_2.out
94: nsize: 2
95: test:
96: suffix: 2
97: args: -mat_type {{aij dense}}
98: test:
99: suffix: 2_scalapack
100: requires: scalapack
101: args: -mat_type scalapack
103: TEST*/