Actual source code: ex162.c
1: static char help[] = "Tests MatShift for SeqAIJ matrices with some missing diagonal entries\n\n";
3: #include <petscmat.h>
5: int main(int argc, char **argv)
6: {
7: Mat A;
8: PetscInt coli[4], row;
9: PetscScalar vali[4];
10: PetscMPIInt size;
12: PetscFunctionBeginUser;
13: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
14: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
15: PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!");
17: PetscCall(MatCreate(PETSC_COMM_SELF, &A));
18: PetscCall(MatSetSizes(A, 4, 4, 4, 4));
19: PetscCall(MatSetType(A, MATSEQAIJ));
20: PetscCall(MatSeqAIJSetPreallocation(A, 4, NULL));
22: row = 0;
23: coli[0] = 1;
24: coli[1] = 3;
25: vali[0] = 1.0;
26: vali[1] = 2.0;
27: PetscCall(MatSetValues(A, 1, &row, 2, coli, vali, ADD_VALUES));
29: row = 1;
30: coli[0] = 0;
31: coli[1] = 1;
32: coli[2] = 2;
33: coli[3] = 3;
34: vali[0] = 3.0;
35: vali[1] = 4.0;
36: vali[2] = 5.0;
37: vali[3] = 6.0;
38: PetscCall(MatSetValues(A, 1, &row, 4, coli, vali, ADD_VALUES));
40: row = 2;
41: coli[0] = 0;
42: coli[1] = 3;
43: vali[0] = 7.0;
44: vali[1] = 8.0;
45: PetscCall(MatSetValues(A, 1, &row, 2, coli, vali, ADD_VALUES));
47: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
48: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
49: PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD));
51: PetscCall(MatShift(A, 0.0));
52: PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD));
54: PetscCall(MatDestroy(&A));
55: PetscCall(PetscFinalize());
56: return 0;
57: }
59: /*TEST
61: test:
63: TEST*/