GCC Code Coverage Report

Directory: ./
File: src/eps/tutorials/ex41.c
Date: 2026-02-22 03:58:10
Exec Total Coverage
Lines: 113 125 90.4%
Functions: 3 3 100.0%
Branches: 334 576 58.0%
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 32 int main(int argc,char **argv)
25 {
26 32 Vec v0,w0; /* initial vectors */
27 32 Mat A; /* operator matrix */
28 32 EPS eps; /* eigenproblem solver context */
29 32 EPSType type;
30 32 PetscInt i,N,m=15,nconv;
31 32 PetscBool twosided;
32 32 PetscReal nrmr,nrml=0.0,re,im,lev;
33 32 PetscScalar *kr,*ki;
34 32 Vec t,*xr,*xi,*yr,*yi;
35 32 PetscMPIInt rank;
36
37
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 32 PetscFunctionBeginUser;
38
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
39
40
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL));
41 32 N = m*(m+1)/2;
42
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
49
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N));
50
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatSetFromOptions(A));
51
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatMarkovModel(m,A));
52
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
59
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSSetOperators(eps,A,NULL));
60
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSSetFromOptions(eps));
67
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatCreateVecs(A,&v0,&w0));
74
14/28
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 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.
 32 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank));
75
1/2
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
 32 if (!rank) {
76
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecSetValue(v0,0,1.0,INSERT_VALUES));
77
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecSetValue(v0,1,1.0,INSERT_VALUES));
78
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecSetValue(v0,2,1.0,INSERT_VALUES));
79
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecSetValue(w0,0,2.0,INSERT_VALUES));
80
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecSetValue(w0,2,0.5,INSERT_VALUES));
81 }
82
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecAssemblyBegin(v0));
83
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecAssemblyBegin(w0));
84
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecAssemblyEnd(v0));
85
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecAssemblyEnd(w0));
86
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSSetInitialSpace(eps,1,&v0));
87
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSSetLeftInitialSpace(eps,1,&w0));
88
89 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
90 Solve the eigensystem
91 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
92
93
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSGetType(eps,&type));
99
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSGetConverged(eps,&nconv));
109
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(PetscMalloc2(nconv,&kr,nconv,&ki));
111
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDuplicateVecs(t,nconv,&xr));
112
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDuplicateVecs(t,nconv,&xi));
113
1/2
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
 32 if (twosided) {
114
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDuplicateVecs(t,nconv,&yr));
115
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDuplicateVecs(t,nconv,&yi));
116 }
117
118
1/2
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
 32 if (nconv>0) {
119 /*
120 Display eigenvalues and relative errors
121 */
122
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(PetscPrintf(PETSC_COMM_WORLD,
123 " k ||Ax-kx|| ||y'A-y'k||\n"
124 " ---------------- ------------------ ------------------\n"));
125
126
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 168 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 136 PetscCall(EPSGetEigenpair(eps,i,&kr[i],&ki[i],xr[i],xi[i]));
132
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 136 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 136 PetscCall(ComputeResidualNorm(A,PETSC_FALSE,kr[i],ki[i],xr[i],xi[i],t,&nrmr));
137
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 136 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 68 re = PetscRealPart(kr[i]);
141 68 im = PetscImaginaryPart(kr[i]);
142 #else
143 68 re = kr[i];
144 68 im = ki[i];
145 #endif
146
6/8
✓ Branch 0 taken 4 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 3 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
 136 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 136 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n"));
150 /*
151 Check bi-orthogonality of eigenvectors
152 */
153
1/2
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
 32 if (twosided) {
154
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecCheckOrthogonality(xr,nconv,yr,nconv,NULL,NULL,&lev));
155
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 32 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(EPSDestroy(&eps));
161
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatDestroy(&A));
162
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDestroy(&v0));
163
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDestroy(&w0));
164
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDestroy(&t));
165
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(PetscFree2(kr,ki));
166
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDestroyVecs(nconv,&xr));
167
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDestroyVecs(nconv,&xi));
168
1/2
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
 32 if (twosided) {
169
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDestroyVecs(nconv,&yr));
170
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(VecDestroyVecs(nconv,&yi));
171 }
172
2/6
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 32 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 32 PetscErrorCode MatMarkovModel(PetscInt m,Mat A)
196 {
197 32 const PetscReal cst = 0.5/(PetscReal)(m-1);
198 32 PetscReal pd,pu;
199 32 PetscInt Istart,Iend,i,j,jmax,ix=0;
200
201
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 32 PetscFunctionBeginUser;
202
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
203
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 512 for (i=1;i<=m;i++) {
204 480 jmax = m-i+1;
205
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 4320 for (j=1;j<=jmax;j++) {
206 3840 ix = ix + 1;
207
2/4
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 8 times.
 3840 if (ix-1<Istart || ix>Iend) continue; /* compute only owned rows */
208
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 3840 if (j!=jmax) {
209 3360 pd = cst*(PetscReal)(i+j-1);
210 /* north */
211
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 3360 if (i==1) PetscCall(MatSetValue(A,ix-1,ix,2*pd,INSERT_VALUES));
212
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 2912 else PetscCall(MatSetValue(A,ix-1,ix,pd,INSERT_VALUES));
213 /* east */
214
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 3360 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 2912 else PetscCall(MatSetValue(A,ix-1,ix+jmax-1,pd,INSERT_VALUES));
216 }
217 /* south */
218 3840 pu = 0.5 - cst*(PetscReal)(i+j-3);
219
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 3840 if (j>1) PetscCall(MatSetValue(A,ix-1,ix-2,pu,INSERT_VALUES));
220 /* west */
221
6/8
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 3840 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 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
225
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 32 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 272 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 136 PetscReal ni,nr;
243 #endif
244
2/2
✓ Branch 0 taken 8 times.
✓ Branch 1 taken 8 times.
 272 PetscErrorCode (*matmult)(Mat,Vec,Vec) = trans? MatMultTranspose: MatMult;
245
246
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 272 PetscFunctionBegin;
247 #if !defined(PETSC_USE_COMPLEX)
248
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
 136 if (ki == 0 || PetscAbsScalar(ki) < PetscAbsScalar(kr*PETSC_MACHINE_EPSILON)) {
249 #endif
250
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 272 PetscCall((*matmult)(A,xr,u));
251
5/8
✓ Branch 0 taken 8 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
 272 if (PetscAbsScalar(kr) > PETSC_MACHINE_EPSILON) PetscCall(VecAXPY(u,-kr,xr));
252
4/6
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
 272 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