xref: /petsc/src/mat/tests/ex203.c (revision 58d68138c660dfb4e9f5b03334792cd4f2ffd7cc)
1 
2 static char help[] = "Tests incorrect use of MatDiagonalSet() for SHELL matrices\n\n";
3 
4 #include <petscmat.h>
5 
6 typedef struct _n_User *User;
7 struct _n_User {
8   Mat B;
9 };
10 
11 static PetscErrorCode MatGetDiagonal_User(Mat A, Vec X) {
12   User user;
13 
14   PetscFunctionBegin;
15   PetscCall(MatShellGetContext(A, &user));
16   PetscCall(MatGetDiagonal(user->B, X));
17   PetscFunctionReturn(0);
18 }
19 
20 int main(int argc, char **args) {
21   const PetscScalar xvals[] = {11, 13};
22   const PetscInt    inds[]  = {0, 1};
23   PetscScalar       avals[] = {2, 3, 5, 7};
24   Mat               A, S;
25   Vec               X, Y;
26   User              user;
27 
28   PetscFunctionBeginUser;
29   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
30   PetscCall(MatCreateSeqAIJ(PETSC_COMM_WORLD, 2, 2, 2, NULL, &A));
31   PetscCall(MatSetUp(A));
32   PetscCall(MatSetValues(A, 2, inds, 2, inds, avals, INSERT_VALUES));
33   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
34   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
35   PetscCall(VecCreateSeq(PETSC_COMM_WORLD, 2, &X));
36   PetscCall(VecSetValues(X, 2, inds, xvals, INSERT_VALUES));
37   PetscCall(VecAssemblyBegin(X));
38   PetscCall(VecAssemblyEnd(X));
39   PetscCall(VecDuplicate(X, &Y));
40 
41   PetscCall(PetscNew(&user));
42   user->B = A;
43 
44   PetscCall(MatCreateShell(PETSC_COMM_WORLD, 2, 2, 2, 2, user, &S));
45   PetscCall(MatShellSetOperation(S, MATOP_GET_DIAGONAL, (void (*)(void))MatGetDiagonal_User));
46   PetscCall(MatSetUp(S));
47 
48   PetscCall(MatShift(S, 42));
49   PetscCall(MatGetDiagonal(S, Y));
50   PetscCall(VecView(Y, PETSC_VIEWER_STDOUT_WORLD));
51   PetscCall(MatDiagonalSet(S, X, ADD_VALUES));
52   PetscCall(MatGetDiagonal(S, Y));
53   PetscCall(VecView(Y, PETSC_VIEWER_STDOUT_WORLD));
54   PetscCall(MatScale(S, 42));
55   PetscCall(MatGetDiagonal(S, Y));
56   PetscCall(VecView(Y, PETSC_VIEWER_STDOUT_WORLD));
57 
58   PetscCall(MatDestroy(&A));
59   PetscCall(MatDestroy(&S));
60   PetscCall(VecDestroy(&X));
61   PetscCall(VecDestroy(&Y));
62   PetscCall(PetscFree(user));
63   PetscCall(PetscFinalize());
64   return 0;
65 }
66 
67 /*TEST
68 
69    test:
70       args: -malloc_dump
71 
72 TEST*/
73