GCC Code Coverage Report


Directory: ./
File: src/eps/tutorials/ex41.c
Date: 2025-10-03 04:28:47
Exec Total Coverage
Lines: 113 125 90.4%
Functions: 3 3 100.0%
Branches: 335 576 58.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[] = "Illustrates the computation of left eigenvectors.\n\n"
12 "The problem is the Markov model as in ex5.c.\n"
13 "The command line options are:\n"
14 " -m <m>, where <m> = number of grid subdivisions in each dimension.\n\n";
15
16 #include <slepceps.h>
17
18 /*
19 User-defined routines
20 */
21 PetscErrorCode MatMarkovModel(PetscInt,Mat);
22 PetscErrorCode ComputeResidualNorm(Mat,PetscBool,PetscScalar,PetscScalar,Vec,Vec,Vec,PetscReal*);
23
24 40 int main(int argc,char **argv)
25 {
26 40 Vec v0,w0; /* initial vectors */
27 40 Mat A; /* operator matrix */
28 40 EPS eps; /* eigenproblem solver context */
29 40 EPSType type;
30 40 PetscInt i,N,m=15,nconv;
31 40 PetscBool twosided;
32 40 PetscReal nrmr,nrml=0.0,re,im,lev;
33 40 PetscScalar *kr,*ki;
34 40 Vec t,*xr,*xi,*yr,*yi;
35 40 PetscMPIInt rank;
36
37
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
40 PetscFunctionBeginUser;
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.
40 PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
39
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.
40 PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL));
41 40 N = m*(m+1)/2;
42
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.
40 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nMarkov Model, N=%" PetscInt_FMT " (m=%" PetscInt_FMT ")\n\n",N,m));
43
44 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
45 Compute the operator matrix that defines the eigensystem, Ax=kx
46 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
47
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.
40 PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
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.
40 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N));
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.
40 PetscCall(MatSetFromOptions(A));
51
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.
40 PetscCall(MatMarkovModel(m,A));
52
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.
40 PetscCall(MatCreateVecs(A,NULL,&t));
53
54 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
55 Create the eigensolver and set various options
56 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
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.
40 PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
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.
40 PetscCall(EPSSetOperators(eps,A,NULL));
60
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.
40 PetscCall(EPSSetProblemType(eps,EPS_NHEP));
61
62 /* use a two-sided algorithm to compute left eigenvectors as well */
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.
40 PetscCall(EPSSetTwoSided(eps,PETSC_TRUE));
64
65 /* allow user to change settings at run time */
66
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.
40 PetscCall(EPSSetFromOptions(eps));
67
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.
40 PetscCall(EPSGetTwoSided(eps,&twosided));
68
69 /*
70 Set the initial vectors. This is optional, if not done the initial
71 vectors are set to random values
72 */
73
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.
40 PetscCall(MatCreateVecs(A,&v0,&w0));
74
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.
40 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank));
75
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
40 if (!rank) {
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.
40 PetscCall(VecSetValue(v0,0,1.0,INSERT_VALUES));
77
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.
40 PetscCall(VecSetValue(v0,1,1.0,INSERT_VALUES));
78
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.
40 PetscCall(VecSetValue(v0,2,1.0,INSERT_VALUES));
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.
40 PetscCall(VecSetValue(w0,0,2.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.
40 PetscCall(VecSetValue(w0,2,0.5,INSERT_VALUES));
81 }
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.
40 PetscCall(VecAssemblyBegin(v0));
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.
40 PetscCall(VecAssemblyBegin(w0));
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.
40 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.
40 PetscCall(VecAssemblyEnd(w0));
86
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.
40 PetscCall(EPSSetInitialSpace(eps,1,&v0));
87
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.
40 PetscCall(EPSSetLeftInitialSpace(eps,1,&w0));
88
89 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
90 Solve the eigensystem
91 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
92
93
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.
40 PetscCall(EPSSolve(eps));
94
95 /*
96 Optional: Get some information from the solver and display it
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.
40 PetscCall(EPSGetType(eps,&type));
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.
40 PetscCall(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type));
100
101 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
102 Display solution and clean up
103 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
104
105 /*
106 Get number of converged approximate eigenpairs
107 */
108
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.
40 PetscCall(EPSGetConverged(eps,&nconv));
109
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.
40 PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of converged eigenpairs: %" PetscInt_FMT "\n\n",nconv));
110
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.
40 PetscCall(PetscMalloc2(nconv,&kr,nconv,&ki));
111
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.
40 PetscCall(VecDuplicateVecs(t,nconv,&xr));
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.
40 PetscCall(VecDuplicateVecs(t,nconv,&xi));
113
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
40 if (twosided) {
114
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.
40 PetscCall(VecDuplicateVecs(t,nconv,&yr));
115
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.
40 PetscCall(VecDuplicateVecs(t,nconv,&yi));
116 }
117
118
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
40 if (nconv>0) {
119 /*
120 Display eigenvalues and relative errors
121 */
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.
40 PetscCall(PetscPrintf(PETSC_COMM_WORLD,
123 " k ||Ax-kx|| ||y'A-y'k||\n"
124 " ---------------- ------------------ ------------------\n"));
125
126
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
210 for (i=0;i<nconv;i++) {
127 /*
128 Get converged eigenpairs: i-th eigenvalue is stored in kr (real part) and
129 ki (imaginary part)
130 */
131
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.
170 PetscCall(EPSGetEigenpair(eps,i,&kr[i],&ki[i],xr[i],xi[i]));
132
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.
170 if (twosided) PetscCall(EPSGetLeftEigenvector(eps,i,yr[i],yi[i]));
133 /*
134 Compute the residual norms associated to each eigenpair
135 */
136
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.
170 PetscCall(ComputeResidualNorm(A,PETSC_FALSE,kr[i],ki[i],xr[i],xi[i],t,&nrmr));
137
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.
170 if (twosided) PetscCall(ComputeResidualNorm(A,PETSC_TRUE,kr[i],ki[i],yr[i],yi[i],t,&nrml));
138
139 #if defined(PETSC_USE_COMPLEX)
140 85 re = PetscRealPart(kr[i]);
141 85 im = PetscImaginaryPart(kr[i]);
142 #else
143 85 re = kr[i];
144 85 im = ki[i];
145 #endif
146
6/8
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
170 if (im!=0.0) PetscCall(PetscPrintf(PETSC_COMM_WORLD," %8f%+8fi %12g %12g\n",(double)re,(double)im,(double)nrmr,(double)nrml));
147
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.
170 else PetscCall(PetscPrintf(PETSC_COMM_WORLD," %12f %12g %12g\n",(double)re,(double)nrmr,(double)nrml));
148 }
149
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.
40 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n"));
150 /*
151 Check bi-orthogonality of eigenvectors
152 */
153
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
40 if (twosided) {
154
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.
40 PetscCall(VecCheckOrthogonality(xr,nconv,yr,nconv,NULL,NULL,&lev));
155
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.
40 if (lev<100*PETSC_MACHINE_EPSILON) PetscCall(PetscPrintf(PETSC_COMM_WORLD," Level of bi-orthogonality of eigenvectors < 100*eps\n\n"));
156 else PetscCall(PetscPrintf(PETSC_COMM_WORLD," Level of bi-orthogonality of eigenvectors: %g\n\n",(double)lev));
157 }
158 }
159
160
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.
40 PetscCall(EPSDestroy(&eps));
161
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.
40 PetscCall(MatDestroy(&A));
162
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.
40 PetscCall(VecDestroy(&v0));
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.
40 PetscCall(VecDestroy(&w0));
164
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.
40 PetscCall(VecDestroy(&t));
165
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.
40 PetscCall(PetscFree2(kr,ki));
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.
40 PetscCall(VecDestroyVecs(nconv,&xr));
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.
40 PetscCall(VecDestroyVecs(nconv,&xi));
168
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
40 if (twosided) {
169
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.
40 PetscCall(VecDestroyVecs(nconv,&yr));
170
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.
40 PetscCall(VecDestroyVecs(nconv,&yi));
171 }
172
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.
40 PetscCall(SlepcFinalize());
173 return 0;
174 }
175
176 /*
177 Matrix generator for a Markov model of a random walk on a triangular grid.
178
179 This subroutine generates a test matrix that models a random walk on a
180 triangular grid. This test example was used by G. W. Stewart ["{SRRIT} - a
181 FORTRAN subroutine to calculate the dominant invariant subspaces of a real
182 matrix", Tech. report. TR-514, University of Maryland (1978).] and in a few
183 papers on eigenvalue problems by Y. Saad [see e.g. LAA, vol. 34, pp. 269-295
184 (1980) ]. These matrices provide reasonably easy test problems for eigenvalue
185 algorithms. The transpose of the matrix is stochastic and so it is known
186 that one is an exact eigenvalue. One seeks the eigenvector of the transpose
187 associated with the eigenvalue unity. The problem is to calculate the steady
188 state probability distribution of the system, which is the eigevector
189 associated with the eigenvalue one and scaled in such a way that the sum all
190 the components is equal to one.
191
192 Note: the code will actually compute the transpose of the stochastic matrix
193 that contains the transition probabilities.
194 */
195 40 PetscErrorCode MatMarkovModel(PetscInt m,Mat A)
196 {
197 40 const PetscReal cst = 0.5/(PetscReal)(m-1);
198 40 PetscReal pd,pu;
199 40 PetscInt Istart,Iend,i,j,jmax,ix=0;
200
201
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
40 PetscFunctionBeginUser;
202
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.
40 PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
203
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
640 for (i=1;i<=m;i++) {
204 600 jmax = m-i+1;
205
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
5400 for (j=1;j<=jmax;j++) {
206 4800 ix = ix + 1;
207
2/4
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
4800 if (ix-1<Istart || ix>Iend) continue; /* compute only owned rows */
208
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
4800 if (j!=jmax) {
209 4200 pd = cst*(PetscReal)(i+j-1);
210 /* north */
211
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.
4200 if (i==1) PetscCall(MatSetValue(A,ix-1,ix,2*pd,INSERT_VALUES));
212
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.
3640 else PetscCall(MatSetValue(A,ix-1,ix,pd,INSERT_VALUES));
213 /* east */
214
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.
4200 if (j==1) PetscCall(MatSetValue(A,ix-1,ix+jmax-1,2*pd,INSERT_VALUES));
215
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.
3640 else PetscCall(MatSetValue(A,ix-1,ix+jmax-1,pd,INSERT_VALUES));
216 }
217 /* south */
218 4800 pu = 0.5 - cst*(PetscReal)(i+j-3);
219
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.
4800 if (j>1) PetscCall(MatSetValue(A,ix-1,ix-2,pu,INSERT_VALUES));
220 /* west */
221
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.
4800 if (i>1) PetscCall(MatSetValue(A,ix-1,ix-jmax-2,pu,INSERT_VALUES));
222 }
223 }
224
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.
40 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
225
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.
40 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
226
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.
8 PetscFunctionReturn(PETSC_SUCCESS);
227 }
228
229 /*
230 ComputeResidualNorm - Computes the norm of the residual vector
231 associated with an eigenpair.
232
233 Input Parameters:
234 trans - whether A' must be used instead of A
235 kr,ki - eigenvalue
236 xr,xi - eigenvector
237 u - work vector
238 */
239 340 PetscErrorCode ComputeResidualNorm(Mat A,PetscBool trans,PetscScalar kr,PetscScalar ki,Vec xr,Vec xi,Vec u,PetscReal *norm)
240 {
241 #if !defined(PETSC_USE_COMPLEX)
242 170 PetscReal ni,nr;
243 #endif
244
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
340 PetscErrorCode (*matmult)(Mat,Vec,Vec) = trans? MatMultTranspose: MatMult;
245
246
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
340 PetscFunctionBegin;
247 #if !defined(PETSC_USE_COMPLEX)
248
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
170 if (ki == 0 || PetscAbsScalar(ki) < PetscAbsScalar(kr*PETSC_MACHINE_EPSILON)) {
249 #endif
250
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.
340 PetscCall((*matmult)(A,xr,u));
251
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.
340 if (PetscAbsScalar(kr) > PETSC_MACHINE_EPSILON) PetscCall(VecAXPY(u,-kr,xr));
252
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.
340 PetscCall(VecNorm(u,NORM_2,norm));
253 #if !defined(PETSC_USE_COMPLEX)
254 } else {
255 PetscCall((*matmult)(A,xr,u));
256 if (SlepcAbsEigenvalue(kr,ki) > PETSC_MACHINE_EPSILON) {
257 PetscCall(VecAXPY(u,-kr,xr));
258 PetscCall(VecAXPY(u,ki,xi));
259 }
260 PetscCall(VecNorm(u,NORM_2,&nr));
261 PetscCall((*matmult)(A,xi,u));
262 if (SlepcAbsEigenvalue(kr,ki) > PETSC_MACHINE_EPSILON) {
263 PetscCall(VecAXPY(u,-kr,xi));
264 PetscCall(VecAXPY(u,-ki,xr));
265 }
266 PetscCall(VecNorm(u,NORM_2,&ni));
267 *norm = SlepcAbsEigenvalue(nr,ni);
268 }
269 #endif
270
6/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 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
68 PetscFunctionReturn(PETSC_SUCCESS);
271 }
272
273 /*TEST
274
275 testset:
276 args: -st_type sinvert -eps_target 1.1 -eps_nev 4
277 filter: grep -v method | sed -e "s/[+-]0\.0*i//g" | sed -e "s/[0-9]\.[0-9]*e[+-]\([0-9]*\)/removed/g"
278 requires: !single
279 output_file: output/ex41_1.out
280 test:
281 suffix: 1
282 args: -eps_type {{power krylovschur}}
283 test:
284 suffix: 1_balance
285 args: -eps_balance {{oneside twoside}} -eps_ncv 17 -eps_krylovschur_locking 0
286 requires: !__float128
287
288 TEST*/
289