GCC Code Coverage Report


Directory: ./
File: src/eps/tutorials/ex5.c
Date: 2025-10-03 04:28:47
Exec Total Coverage
Lines: 70 74 94.6%
Functions: 2 2 100.0%
Branches: 207 338 61.2%

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[] = "Eigenvalue problem associated with a Markov model of a random walk on a triangular grid.\n\n"
12 "It is a standard nonsymmetric eigenproblem with real eigenvalues and the rightmost eigenvalue is known to be 1.\n"
13 "This example illustrates how the user can set the initial vector.\n\n"
14 "The command line options are:\n"
15 " -m <m>, where <m> = number of grid subdivisions in each dimension.\n\n";
16
17 #include <slepceps.h>
18
19 /*
20 User-defined routines
21 */
22 PetscErrorCode MatMarkovModel(PetscInt m,Mat A);
23
24 90 int main(int argc,char **argv)
25 {
26 90 Vec v0; /* initial vector */
27 90 Mat A; /* operator matrix */
28 90 EPS eps; /* eigenproblem solver context */
29 90 EPSType type;
30 90 EPSStop stop;
31 90 PetscReal thres;
32 90 PetscInt N,m=15,nev;
33 90 PetscMPIInt rank;
34 90 PetscBool terse;
35
36
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
90 PetscFunctionBeginUser;
37
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.
90 PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
38
39
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.
90 PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL));
40 90 N = m*(m+1)/2;
41
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.
90 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nMarkov Model, N=%" PetscInt_FMT " (m=%" PetscInt_FMT ")\n\n",N,m));
42
43 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
44 Compute the operator matrix that defines the eigensystem, Ax=kx
45 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
46
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.
90 PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
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.
90 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N));
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.
90 PetscCall(MatSetFromOptions(A));
50
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.
90 PetscCall(MatMarkovModel(m,A));
51
52 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
53 Create the eigensolver and set various options
54 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
55
56 /*
57 Create eigensolver context
58 */
59
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.
90 PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
60
61 /*
62 Set operators. In this case, it is a standard eigenvalue problem
63 */
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.
90 PetscCall(EPSSetOperators(eps,A,NULL));
65
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.
90 PetscCall(EPSSetProblemType(eps,EPS_NHEP));
66
67 /*
68 Set solver parameters at runtime
69 */
70
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.
90 PetscCall(EPSSetFromOptions(eps));
71
72 /*
73 Set the initial vector. This is optional, if not done the initial
74 vector is set to random values
75 */
76
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.
90 PetscCall(MatCreateVecs(A,&v0,NULL));
77
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.
90 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank));
78
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
90 if (!rank) {
79
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.
50 PetscCall(VecSetValue(v0,0,1.0,INSERT_VALUES));
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.
50 PetscCall(VecSetValue(v0,1,1.0,INSERT_VALUES));
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.
50 PetscCall(VecSetValue(v0,2,1.0,INSERT_VALUES));
82 }
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.
90 PetscCall(VecAssemblyBegin(v0));
84
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.
90 PetscCall(VecAssemblyEnd(v0));
85
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.
90 PetscCall(EPSSetInitialSpace(eps,1,&v0));
86
87 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
88 Solve the eigensystem
89 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
90
91
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.
90 PetscCall(EPSSolve(eps));
92
93 /*
94 Optional: Get some information from the solver and display it
95 */
96
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.
90 PetscCall(EPSGetType(eps,&type));
97
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.
90 PetscCall(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type));
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.
90 PetscCall(EPSGetStoppingTest(eps,&stop));
99
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
90 if (stop!=EPS_STOP_THRESHOLD) {
100
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.
80 PetscCall(EPSGetDimensions(eps,&nev,NULL,NULL));
101
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.
80 PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of requested eigenvalues: %" PetscInt_FMT "\n",nev));
102 } else {
103
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(EPSGetThreshold(eps,&thres,NULL));
104
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(PetscPrintf(PETSC_COMM_WORLD," Using threshold: %.4g\n",(double)thres));
105 }
106
107 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
108 Display solution and clean up
109 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
110
111 /* show detailed info unless -terse option is given by user */
112
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.
90 PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
113
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.
90 if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
114 else {
115 PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL));
116 PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD));
117 PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD));
118 PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
119 }
120
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.
90 PetscCall(EPSDestroy(&eps));
121
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.
90 PetscCall(MatDestroy(&A));
122
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.
90 PetscCall(VecDestroy(&v0));
123
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.
90 PetscCall(SlepcFinalize());
124 return 0;
125 }
126
127 /*
128 Matrix generator for a Markov model of a random walk on a triangular grid.
129
130 This subroutine generates a test matrix that models a random walk on a
131 triangular grid. This test example was used by G. W. Stewart ["{SRRIT} - a
132 FORTRAN subroutine to calculate the dominant invariant subspaces of a real
133 matrix", Tech. report. TR-514, University of Maryland (1978).] and in a few
134 papers on eigenvalue problems by Y. Saad [see e.g. LAA, vol. 34, pp. 269-295
135 (1980) ]. These matrices provide reasonably easy test problems for eigenvalue
136 algorithms. The transpose of the matrix is stochastic and so it is known
137 that one is an exact eigenvalue. One seeks the eigenvector of the transpose
138 associated with the eigenvalue unity. The problem is to calculate the steady
139 state probability distribution of the system, which is the eigevector
140 associated with the eigenvalue one and scaled in such a way that the sum all
141 the components is equal to one.
142
143 Note: the code will actually compute the transpose of the stochastic matrix
144 that contains the transition probabilities.
145 */
146 90 PetscErrorCode MatMarkovModel(PetscInt m,Mat A)
147 {
148 90 const PetscReal cst = 0.5/(PetscReal)(m-1);
149 90 PetscReal pd,pu;
150 90 PetscInt Istart,Iend,i,j,jmax,ix=0;
151
152
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
90 PetscFunctionBeginUser;
153
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.
90 PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
154
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
1440 for (i=1;i<=m;i++) {
155 1350 jmax = m-i+1;
156
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
12150 for (j=1;j<=jmax;j++) {
157 10800 ix = ix + 1;
158
4/4
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
10800 if (ix-1<Istart || ix>Iend) continue; /* compute only owned rows */
159
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
6000 if (j!=jmax) {
160 5250 pd = cst*(PetscReal)(i+j-1);
161 /* north */
162
6/8
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ 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.
5250 if (i==1) PetscCall(MatSetValue(A,ix-1,ix,2*pd,INSERT_VALUES));
163
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.
4550 else PetscCall(MatSetValue(A,ix-1,ix,pd,INSERT_VALUES));
164 /* east */
165
6/8
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ 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.
5250 if (j==1) PetscCall(MatSetValue(A,ix-1,ix+jmax-1,2*pd,INSERT_VALUES));
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.
4550 else PetscCall(MatSetValue(A,ix-1,ix+jmax-1,pd,INSERT_VALUES));
167 }
168 /* south */
169 6000 pu = 0.5 - cst*(PetscReal)(i+j-3);
170
6/8
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ 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.
6000 if (j>1) PetscCall(MatSetValue(A,ix-1,ix-2,pu,INSERT_VALUES));
171 /* west */
172
6/8
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ 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.
10800 if (i>1) PetscCall(MatSetValue(A,ix-1,ix-jmax-2,pu,INSERT_VALUES));
173 }
174 }
175
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.
90 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
176
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.
90 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
177
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.
18 PetscFunctionReturn(PETSC_SUCCESS);
178 }
179
180 /*TEST
181
182 test:
183 suffix: 1
184 nsize: 2
185 args: -eps_largest_real -eps_nev 4 -eps_two_sided {{0 1}} -eps_krylovschur_locking {{0 1}} -ds_parallel synchronized -terse
186 filter: sed -e "s/90424/90423/" | sed -e "s/85715/85714/"
187
188 test:
189 suffix: 2
190 args: -eps_threshold_relative .9 -eps_ncv 10 -terse
191 filter: sed -e "s/-0.85714/0.85714/" | sed -e "s/90424/90423/" | sed -e "s/-1.00000, 1.00000/1.00000, -1.00000/" | sed -e "s/-0.97137, 0.97137/0.97137, -0.97137/" | sed -e "s/-0.90423, 0.90423/0.90423, -0.90423/"
192
193 TEST*/
194