xref: /petsc/src/mat/tests/ex3.c (revision ebead697dbf761eb322f829370bbe90b3bd93fa3)
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.; v[1] = 2.; v[2] = -1.;
31   for (i=1; i<n-1; i++) {
32     midx[0] = i-1; midx[1] = i; midx[2] = i+1;
33     PetscCall(MatSetValues(C,1,&i,3,midx,v,INSERT_VALUES));
34   }
35   i    = 0; midx[0] = 0; midx[1] = 1;
36   v[0] = 2.0; v[1] = -1.;
37   PetscCall(MatSetValues(C,1,&i,2,midx,v,INSERT_VALUES));
38   i    = n-1; midx[0] = n-2; midx[1] = n-1;
39   v[0] = -1.0; v[1] = 2.;
40   PetscCall(MatSetValues(C,1,&i,2,midx,v,INSERT_VALUES));
41 
42   PetscCall(MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY));
43   PetscCall(MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY));
44 
45   PetscCall(MatMult(C,u,b));
46 
47   for (i=0; i<n; i++) {
48     PetscCall(MatSOR(C,b,omega,SOR_FORWARD_SWEEP,0.0,1,1,x));
49     PetscCall(VecWAXPY(e,-1.0,x,u));
50     PetscCall(VecNorm(e,NORM_2,&norm));
51     PetscCall(PetscPrintf(PETSC_COMM_SELF,"2-norm of error %g\n",(double)norm));
52   }
53   PetscCall(MatDestroy(&C));
54   PetscCall(VecDestroy(&x));
55   PetscCall(VecDestroy(&b));
56   PetscCall(VecDestroy(&u));
57   PetscCall(VecDestroy(&e));
58   PetscCall(PetscFinalize());
59   return 0;
60 }
61 
62 /*TEST
63 
64    test:
65 
66 TEST*/
67