xref: /petsc/src/mat/tests/ex3.c (revision df4cd43f92eaa320656440c40edb1046daee8f75)
1 
2 static char help[] = "Tests relaxation for dense matrices.\n\n";
3 
4 #include <petscmat.h>
5 
6 int main(int argc, char **args)
7 {
8   Mat         C;
9   Vec         u, x, b, e;
10   PetscInt    i, n = 10, midx[3];
11   PetscScalar v[3];
12   PetscReal   omega = 1.0, norm;
13 
14   PetscFunctionBeginUser;
15   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
16   PetscCall(PetscOptionsGetReal(NULL, NULL, "-omega", &omega, NULL));
17   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
18 
19   PetscCall(MatCreate(PETSC_COMM_SELF, &C));
20   PetscCall(MatSetSizes(C, n, n, n, n));
21   PetscCall(MatSetType(C, MATSEQDENSE));
22   PetscCall(MatSetUp(C));
23   PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &b));
24   PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &x));
25   PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &u));
26   PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &e));
27   PetscCall(VecSet(u, 1.0));
28   PetscCall(VecSet(x, 0.0));
29 
30   v[0] = -1.;
31   v[1] = 2.;
32   v[2] = -1.;
33   for (i = 1; i < n - 1; i++) {
34     midx[0] = i - 1;
35     midx[1] = i;
36     midx[2] = i + 1;
37     PetscCall(MatSetValues(C, 1, &i, 3, midx, v, INSERT_VALUES));
38   }
39   i       = 0;
40   midx[0] = 0;
41   midx[1] = 1;
42   v[0]    = 2.0;
43   v[1]    = -1.;
44   PetscCall(MatSetValues(C, 1, &i, 2, midx, v, INSERT_VALUES));
45   i       = n - 1;
46   midx[0] = n - 2;
47   midx[1] = n - 1;
48   v[0]    = -1.0;
49   v[1]    = 2.;
50   PetscCall(MatSetValues(C, 1, &i, 2, midx, v, INSERT_VALUES));
51 
52   PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY));
53   PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY));
54 
55   PetscCall(MatMult(C, u, b));
56 
57   for (i = 0; i < n; i++) {
58     PetscCall(MatSOR(C, b, omega, SOR_FORWARD_SWEEP, 0.0, 1, 1, x));
59     PetscCall(VecWAXPY(e, -1.0, x, u));
60     PetscCall(VecNorm(e, NORM_2, &norm));
61     PetscCall(PetscPrintf(PETSC_COMM_SELF, "2-norm of error %g\n", (double)norm));
62   }
63   PetscCall(MatDestroy(&C));
64   PetscCall(VecDestroy(&x));
65   PetscCall(VecDestroy(&b));
66   PetscCall(VecDestroy(&u));
67   PetscCall(VecDestroy(&e));
68   PetscCall(PetscFinalize());
69   return 0;
70 }
71 
72 /*TEST
73 
74    test:
75 
76 TEST*/
77