GCC Code Coverage Report

Directory: ./
File: src/eps/tutorials/ex57.c
Date: 2025-10-04 04:19:13
Exec Total Coverage
Lines: 82 87 94.3%
Functions: 2 2 100.0%
Branches: 240 366 65.6%
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[] = "Another eigenvalue problem with Hamiltonian structure.\n\n"
12 "Position and velocity control of a string of high-speed vehicles.\n"
13 "The command line options are:\n"
14 " -m <m>, where <n> = number of vehicles.\n\n";
15
16 #include <slepceps.h>
17
18 /*
19 This example computes eigenvalues of a matrix
20
21 H = [ A N
22 K -A^* ],
23
24 where A, N and K are n-by-n matrices, with n=2*m-1 and
25
26 N = diag( 1, 0, 1, 0,..., 0, 1)
27 K = diag( 0,10, 0,10,...,10, 0)
28 A = tridiag(b,d,c)
29 d = [-1, 0,-1, 0,...,-1, 0,-1]
30 b = [ 1, 0, 1, 0,..., 1, 0]
31 c = [ 0,-1, 0,-1,..., 0,-1]
32
33 References:
34
35 [1] W.R. Ferng, Wen-Wei Lin, Chern-Shuh Wang, The shift-inverted J-Lanczos algorithm
36 for the numerical solutions of large sparse algebraic Riccati equations,
37 Computers & Mathematics with Applications, Volume 33, Issue 10, 1997,
38 https://doi.org/10.1016/S0898-1221(97)00074-6.
39
40 */
41
42 48360 static PetscErrorCode eigenCompare (PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
43 {
44 48360 PetscReal abs1, abs2, tol=1e-12, r_ar, r_ai, r_br, r_bi;
45
46
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 48360 PetscFunctionBegin;
47 #if defined(PETSC_USE_COMPLEX)
48 37840 r_ar = PetscRealPart(ar);
49 37840 r_ai = PetscImaginaryPart(ar);
50 37840 r_br = PetscRealPart(br);
51 37840 r_bi = PetscImaginaryPart(br);
52 #else
53 10520 r_ar = ar;
54 10520 r_ai = ai;
55 10520 r_br = br;
56 10520 r_bi = bi;
57 #endif
58
32/32
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
✓ Branch 4 taken 10 times.
✓ Branch 5 taken 10 times.
✓ Branch 6 taken 10 times.
✓ Branch 7 taken 10 times.
✓ Branch 8 taken 9 times.
✓ Branch 9 taken 10 times.
✓ Branch 10 taken 9 times.
✓ Branch 11 taken 10 times.
✓ Branch 12 taken 10 times.
✓ Branch 13 taken 10 times.
✓ Branch 14 taken 10 times.
✓ Branch 15 taken 10 times.
✓ Branch 16 taken 3 times.
✓ Branch 17 taken 7 times.
✓ Branch 18 taken 3 times.
✓ Branch 19 taken 7 times.
✓ Branch 20 taken 7 times.
✓ Branch 21 taken 7 times.
✓ Branch 22 taken 3 times.
✓ Branch 23 taken 7 times.
✓ Branch 24 taken 3 times.
✓ Branch 25 taken 7 times.
✓ Branch 26 taken 3 times.
✓ Branch 27 taken 7 times.
✓ Branch 28 taken 3 times.
✓ Branch 29 taken 7 times.
✓ Branch 30 taken 7 times.
✓ Branch 31 taken 7 times.
 78434 if (PetscAbs(PetscAbs(r_ar)-PetscAbs(r_br)) < tol && PetscAbs(PetscAbs(r_ai)-PetscAbs(r_bi)) < tol) {
59 /* Negative real part first*/
60
4/4
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 5 times.
 660 if (r_ar<0 && r_br>0)
61 334 *res = -1;
62
4/4
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 5 times.
 326 else if (r_ar>0 && r_br<0)
63 86 *res = 1;
64 /* Positive imaginary part first*/
65
3/4
✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 5 times.
✗ Branch 3 not taken.
 240 else if (r_ai>0 && r_bi<0)
66 205 *res = -1;
67
2/4
✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5 times.
✗ Branch 3 not taken.
 35 else if (r_ai<0 && r_bi>0)
68 35 *res = 1;
69 }
70 else {
71 47700 abs1 = SlepcAbs(r_ar,r_ai);
72 47700 abs2 = SlepcAbs(r_br,r_bi);
73
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
 47700 if (abs1 > abs2)
74 41359 *res = -1;
75
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
 6341 else if (abs1 < abs2)
76 6341 *res = 1;
77 else
78 ✗ *res = 0;
79 }
80
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.
 48360 PetscFunctionReturn(PETSC_SUCCESS);
81 }
82
83 50 int main(int argc,char **argv)
84 {
85 50 Mat H,A,N,K; /* problem matrices */
86 50 EPS eps; /* eigenproblem solver context */
87 50 PetscInt m=12,n,Istart,Iend,i;
88 50 PetscBool terse, sort_hamilt;
89
90
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
 50 PetscFunctionBeginUser;
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.
 50 PetscCall(SlepcInitialize(&argc,&argv,(char*)0,help));
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.
 50 PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL));
94 50 n = 2*m-1;
95
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(PetscPrintf(PETSC_COMM_WORLD,"\nHamiltonian eigenproblem, n=%" PetscInt_FMT
96 " (m=%" PetscInt_FMT " vehicles)\n\n",n,m));
97
98 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
99 Compute the problem matrices A, N and K
100 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
101
102
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(MatCreate(PETSC_COMM_WORLD,&A));
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.
 50 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n));
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.
 50 PetscCall(MatSetFromOptions(A));
105
106
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(MatCreate(PETSC_COMM_WORLD,&N));
107
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(MatSetSizes(N,PETSC_DECIDE,PETSC_DECIDE,n,n));
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.
 50 PetscCall(MatSetFromOptions(N));
109
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.
 50 PetscCall(MatCreate(PETSC_COMM_WORLD,&K));
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.
 50 PetscCall(MatSetSizes(K,PETSC_DECIDE,PETSC_DECIDE,n,n));
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.
 50 PetscCall(MatSetFromOptions(K));
113
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.
 50 PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
115
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
 855 for (i=Istart;i<Iend;i++) {
116
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
 805 if (i%2==0)
117
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.
 420 PetscCall(MatSetValue(A,i,i,-1.0,INSERT_VALUES));
118 else {
119
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.
 385 if (i>0) PetscCall(MatSetValue(A,i,i-1,1.0,INSERT_VALUES));
120
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.
 805 if (i<n-1) PetscCall(MatSetValue(A,i,i+1,-1.0,INSERT_VALUES));
121 }
122 }
123
124
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(MatGetOwnershipRange(N,&Istart,&Iend));
125
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
 470 for (i=Istart+Istart%2;i<Iend;i+=2) {
126
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.
 420 PetscCall(MatSetValue(N,i,i,1.0,INSERT_VALUES));
127 }
128
129
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(MatGetOwnershipRange(K,&Istart,&Iend));
130
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
 435 for (i=Istart+(Istart+1)%2;i<Iend;i+=2) {
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.
 385 PetscCall(MatSetValue(K,i,i,10.0,INSERT_VALUES));
132 }
133
134
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(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
135
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(MatAssemblyBegin(N,MAT_FINAL_ASSEMBLY));
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.
 50 PetscCall(MatAssemblyBegin(K,MAT_FINAL_ASSEMBLY));
137
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(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
138
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(MatAssemblyEnd(N,MAT_FINAL_ASSEMBLY));
139
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(MatAssemblyEnd(K,MAT_FINAL_ASSEMBLY));
140
141
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(MatCreateHamiltonian(A,N,K,&H));
142
143 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
144 Create the eigensolver and set various options
145 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
146
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.
 50 PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
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.
 50 PetscCall(EPSSetOperators(eps,H,NULL));
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.
 50 PetscCall(EPSSetProblemType(eps,EPS_HAMILT));
150
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.
 50 PetscCall(PetscOptionsHasName(NULL,NULL,"-sort_hamilt",&sort_hamilt));
152
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 5 times.
 50 if (sort_hamilt) {
153 /* Adjust ordering of non-hermitian solver so that it is the same as in as in EPS_HAMILT solver */
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.
 30 PetscCall(EPSSetEigenvalueComparison(eps,eigenCompare,NULL));
155
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.
 30 PetscCall(EPSSetWhichEigenpairs(eps,EPS_WHICH_USER));
156 }
157
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(EPSSetFromOptions(eps));
158
159 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
160 Solve the eigensystem and display solution
161 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
162
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.
 50 PetscCall(EPSSolve(eps));
164
165 /* show detailed info unless -terse option is given by user */
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.
 50 PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
167
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.
 50 if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
168 else {
169 ✗ PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL));
170 ✗ PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD));
171 ✗ PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD));
172 ✗ PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
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.
 50 PetscCall(EPSDestroy(&eps));
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.
 50 PetscCall(MatDestroy(&A));
177
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(MatDestroy(&N));
178
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(MatDestroy(&K));
179
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(MatDestroy(&H));
180
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.
 50 PetscCall(SlepcFinalize());
181 return 0;
182 }
183
184 /*TEST
185
186 testset:
187 args: -eps_nev 8 -eps_ncv 28 -terse
188 nsize: {{1 2}}
189 output_file: output/ex57_1.out
190 test:
191 requires: double !complex
192 suffix: 1
193 test:
194 requires: double !complex
195 args: -eps_non_hermitian -sort_hamilt
196 suffix: 1_nhep
197 test:
198 requires: double complex
199 TODO: no support for complex scalars yet
200 suffix: 1_complex
201 test:
202 requires: double complex
203 args: -eps_non_hermitian -sort_hamilt
204 suffix: 1_complex_nhep
205
206 testset:
207 args: -eps_nev 4 -eps_smallest_magnitude -eps_ncv 32 -terse
208 output_file: output/ex57_2.out
209 test:
210 requires: double !complex
211 suffix: 2
212
213 TEST*/
214