Line | Branch | Exec | Source |
---|---|---|---|
1 | /* | ||
2 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
3 | SLEPc - Scalable Library for Eigenvalue Problem Computations | ||
4 | Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain | ||
5 | |||
6 | This file is part of SLEPc. | ||
7 | SLEPc is distributed under a 2-clause BSD license (see LICENSE). | ||
8 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
9 | */ | ||
10 | |||
11 | static char help[] = "Solves the same eigenproblem as in example ex2, but using a shell matrix.\n\n" | ||
12 | "The problem is a standard symmetric eigenproblem corresponding to the 2-D Laplacian operator.\n\n" | ||
13 | "The command line options are:\n" | ||
14 | " -n <n>, where <n> = number of grid subdivisions in both x and y dimensions.\n\n"; | ||
15 | |||
16 | #include <slepceps.h> | ||
17 | |||
18 | /* | ||
19 | User-defined routines | ||
20 | */ | ||
21 | PetscErrorCode MatMult_Laplacian2D(Mat A,Vec x,Vec y); | ||
22 | PetscErrorCode MatGetDiagonal_Laplacian2D(Mat A,Vec diag); | ||
23 | |||
24 | 20 | int main(int argc,char **argv) | |
25 | { | ||
26 | 20 | Mat A; /* operator matrix */ | |
27 | 20 | EPS eps; /* eigenproblem solver context */ | |
28 | 20 | EPSType type; | |
29 | 20 | PetscMPIInt size; | |
30 | 20 | PetscInt N,n=10,nev; | |
31 | 20 | PetscBool terse; | |
32 | |||
33 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
20 | PetscFunctionBeginUser; |
34 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
35 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
20 | PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); |
36 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
20 | PetscCheck(size==1,PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This is a uniprocessor example only"); |
37 | |||
38 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); |
39 | 20 | N = n*n; | |
40 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n2-D Laplacian Eigenproblem (matrix-free version), N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid)\n\n",N,n,n)); |
41 | |||
42 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
43 | Create the operator matrix that defines the eigensystem, Ax=kx | ||
44 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
45 | |||
46 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(MatCreateShell(PETSC_COMM_WORLD,N,N,N,N,&n,&A)); |
47 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(MatShellSetOperation(A,MATOP_MULT,(PetscErrorCodeFn*)MatMult_Laplacian2D)); |
48 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(MatShellSetOperation(A,MATOP_MULT_TRANSPOSE,(PetscErrorCodeFn*)MatMult_Laplacian2D)); |
49 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(MatShellSetOperation(A,MATOP_GET_DIAGONAL,(PetscErrorCodeFn*)MatGetDiagonal_Laplacian2D)); |
50 | |||
51 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
52 | Create the eigensolver and set various options | ||
53 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
54 | |||
55 | /* | ||
56 | Create eigensolver context | ||
57 | */ | ||
58 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps)); |
59 | |||
60 | /* | ||
61 | Set operators. In this case, it is a standard eigenvalue problem | ||
62 | */ | ||
63 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSSetOperators(eps,A,NULL)); |
64 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSSetProblemType(eps,EPS_HEP)); |
65 | |||
66 | /* | ||
67 | Set solver parameters at runtime | ||
68 | */ | ||
69 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSSetFromOptions(eps)); |
70 | |||
71 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
72 | Solve the eigensystem | ||
73 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
74 | |||
75 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSSolve(eps)); |
76 | |||
77 | /* | ||
78 | Optional: Get some information from the solver and display it | ||
79 | */ | ||
80 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSGetType(eps,&type)); |
81 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type)); |
82 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSGetDimensions(eps,&nev,NULL,NULL)); |
83 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of requested eigenvalues: %" PetscInt_FMT "\n",nev)); |
84 | |||
85 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
86 | Display solution and clean up | ||
87 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
88 | |||
89 | /* show detailed info unless -terse option is given by user */ | ||
90 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
91 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
20 | if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL)); |
92 | else { | ||
93 | ✗ | PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL)); | |
94 | ✗ | PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD)); | |
95 | ✗ | PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD)); | |
96 | ✗ | PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); | |
97 | } | ||
98 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(EPSDestroy(&eps)); |
99 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(MatDestroy(&A)); |
100 |
3/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
20 | PetscCall(SlepcFinalize()); |
101 | return 0; | ||
102 | } | ||
103 | |||
104 | /* | ||
105 | Compute the matrix vector multiplication y<---T*x where T is a nx by nx | ||
106 | tridiagonal matrix with DD on the diagonal, DL on the subdiagonal, and | ||
107 | DU on the superdiagonal. | ||
108 | */ | ||
109 | 299160 | static void tv(int nx,const PetscScalar *x,PetscScalar *y) | |
110 | { | ||
111 | 299160 | PetscScalar dd,dl,du; | |
112 | 299160 | int j; | |
113 | |||
114 | 299160 | dd = 4.0; | |
115 | 299160 | dl = -1.0; | |
116 | 299160 | du = -1.0; | |
117 | |||
118 | 299160 | y[0] = dd*x[0] + du*x[1]; | |
119 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
17870040 | for (j=1;j<nx-1;j++) |
120 | 17570880 | y[j] = dl*x[j-1] + dd*x[j] + du*x[j+1]; | |
121 | 299160 | y[nx-1] = dl*x[nx-2] + dd*x[nx-1]; | |
122 | 299160 | } | |
123 | |||
124 | /* | ||
125 | Matrix-vector product subroutine for the 2D Laplacian. | ||
126 | |||
127 | The matrix used is the 2 dimensional discrete Laplacian on unit square with | ||
128 | zero Dirichlet boundary condition. | ||
129 | |||
130 | Computes y <-- A*x, where A is the block tridiagonal matrix | ||
131 | |||
132 | | T -I | | ||
133 | |-I T -I | | ||
134 | A = | -I T | | ||
135 | | ... -I| | ||
136 | | -I T| | ||
137 | |||
138 | The subroutine TV is called to compute y<--T*x. | ||
139 | */ | ||
140 | 8836 | PetscErrorCode MatMult_Laplacian2D(Mat A,Vec x,Vec y) | |
141 | { | ||
142 | 8836 | void *ctx; | |
143 | 8836 | int nx,lo,i,j; | |
144 | 8836 | const PetscScalar *px; | |
145 | 8836 | PetscScalar *py; | |
146 | |||
147 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
8836 | PetscFunctionBeginUser; |
148 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
8836 | PetscCall(MatShellGetContext(A,&ctx)); |
149 | 8836 | nx = *(int*)ctx; | |
150 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
8836 | PetscCall(VecGetArrayRead(x,&px)); |
151 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
8836 | PetscCall(VecGetArray(y,&py)); |
152 | |||
153 | 8836 | tv(nx,&px[0],&py[0]); | |
154 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
316832 | for (i=0;i<nx;i++) py[i] -= px[nx+i]; |
155 | |||
156 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
290324 | for (j=2;j<nx;j++) { |
157 | 281488 | lo = (j-1)*nx; | |
158 | 281488 | tv(nx,&px[lo],&py[lo]); | |
159 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
18133856 | for (i=0;i<nx;i++) py[lo+i] -= px[lo-nx+i] + px[lo+nx+i]; |
160 | } | ||
161 | |||
162 | 8836 | lo = (nx-1)*nx; | |
163 | 8836 | tv(nx,&px[lo],&py[lo]); | |
164 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
316832 | for (i=0;i<nx;i++) py[lo+i] -= px[lo-nx+i]; |
165 | |||
166 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
8836 | PetscCall(VecRestoreArrayRead(x,&px)); |
167 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
8836 | PetscCall(VecRestoreArray(y,&py)); |
168 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
1768 | PetscFunctionReturn(PETSC_SUCCESS); |
169 | } | ||
170 | |||
171 | 10 | PetscErrorCode MatGetDiagonal_Laplacian2D(Mat A,Vec diag) | |
172 | { | ||
173 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBeginUser; |
174 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10 | PetscCall(VecSet(diag,4.0)); |
175 |
5/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
2 | PetscFunctionReturn(PETSC_SUCCESS); |
176 | } | ||
177 | |||
178 | /*TEST | ||
179 | |||
180 | test: | ||
181 | suffix: 1 | ||
182 | args: -n 72 -eps_nev 4 -eps_ncv 20 -terse | ||
183 | requires: !single | ||
184 | |||
185 | test: | ||
186 | suffix: 2 | ||
187 | args: -eps_nev 4 -st_type sinvert -st_ksp_type cg -st_pc_type jacobi -terse | ||
188 | requires: !single | ||
189 | |||
190 | TEST*/ | ||
191 |