GCC Code Coverage Report


Directory: ./
File: src/sys/classes/st/impls/filter/chebyshev.c
Date: 2025-10-04 04:19:13
Exec Total Coverage
Lines: 169 173 97.7%
Functions: 9 9 100.0%
Branches: 292 471 62.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 A polynomial filter based on a truncated Chebyshev series.
12 */
13
14 #include <slepc/private/stimpl.h>
15 #include "filter.h"
16
17 #define DAMPING_LANCZOS_POW 2
18
19 /*
20 Gateway to Chebyshev for evaluating y=p(A)*x
21 (Clenshaw with a vector)
22 */
23 6681 static PetscErrorCode MatMult_Chebyshev(Mat A,Vec x,Vec y)
24 {
25 6681 ST st;
26 6681 ST_FILTER *ctx;
27 6681 CHEBYSHEV_CTX cheby;
28 6681 PetscInt degree,i,w0=0,w1=1;
29 6681 PetscReal s1,s2;
30
31
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
6681 PetscFunctionBegin;
32
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.
6681 PetscCall(MatShellGetContext(A,&st));
33 6681 ctx = (ST_FILTER*)st->data;
34 6681 cheby = (CHEBYSHEV_CTX)ctx->data;
35 6681 degree = ctx->polyDegree;
36 6681 s1 = 4/(ctx->right - ctx->left);
37 6681 s2 = -2*(ctx->right+ctx->left)/(ctx->right-ctx->left);
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.
6681 PetscCall(VecSet(st->work[w0],0.0));
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.
6681 PetscCall(VecCopy(x,y));
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.
6681 PetscCall(VecScale(y,cheby->damping_coeffs[degree]*cheby->coeffs[degree]));
42
43
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
832881 for (i=0;i<degree;i++) {
44 /* swap work vectors, save a copy of current y */
45 826200 w0 = 1-w0;
46 826200 w1 = 1-w1;
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.
826200 PetscCall(VecCopy(y,st->work[w0]));
48 /* compute next y */
49
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
826200 if (i == degree-1) {
50 6681 s1 = s1/2;
51 6681 s2 = s2/2;
52 }
53
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.
826200 PetscCall(VecAXPBYPCZ(y,cheby->damping_coeffs[degree-i-1]*cheby->coeffs[degree-i-1],-1,s2,x,st->work[w1]));
54
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.
826200 PetscCall(MatMult(ctx->T,st->work[w0],st->work[w1]));
55
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.
826200 PetscCall(VecAXPY(y,s1,st->work[w1]));
56 }
57
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.
1263 PetscFunctionReturn(PETSC_SUCCESS);
58 }
59
60 /*
61 Gateway to Chebyshev for evaluating C=p(A)*B
62 (Clenshaw with a matrix)
63 */
64 512 static PetscErrorCode MatMatMult_Chebyshev(Mat A,Mat B,Mat C,void *pctx)
65 {
66 512 ST st;
67 512 ST_FILTER *ctx;
68 512 CHEBYSHEV_CTX cheby;
69 512 PetscInt degree,i,m1,m2,w0=0,w1=1;
70 512 PetscReal s1,s2;
71
72
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
512 PetscFunctionBegin;
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.
512 PetscCall(MatShellGetContext(A,&st));
74 512 ctx = (ST_FILTER*)st->data;
75 512 cheby = (CHEBYSHEV_CTX)ctx->data;
76 512 degree = ctx->polyDegree;
77 512 s1 = 4/(ctx->right - ctx->left);
78 512 s2 = -2*(ctx->right+ctx->left)/(ctx->right-ctx->left);
79
80
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
512 if (ctx->nW) { /* check if work matrices must be resized */
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.
498 PetscCall(MatGetSize(B,NULL,&m1));
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.
498 PetscCall(MatGetSize(ctx->W[0],NULL,&m2));
83
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
498 if (m1!=m2) {
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.
57 PetscCall(MatDestroyMatrices(ctx->nW,&ctx->W));
85 57 ctx->nW = 0;
86 }
87 }
88
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
512 if (!ctx->nW) { /* allocate work matrices */
89 71 ctx->nW = 2;
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.
71 PetscCall(PetscMalloc1(ctx->nW,&ctx->W));
91
7/8
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 2 times.
213 for (i=0;i<ctx->nW;i++) PetscCall(MatDuplicate(B,MAT_DO_NOT_COPY_VALUES,&ctx->W[i]));
92 }
93
94
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.
512 PetscCall(MatScale(ctx->W[w0],0));
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.
512 PetscCall(MatCopy(B,C,SAME_NONZERO_PATTERN));
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.
512 PetscCall(MatScale(C,cheby->damping_coeffs[degree]*cheby->coeffs[degree]));
97
98
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
102912 for (i=0;i<degree;i++) {
99 /* swap work matrices, save a copy of current C */
100 102400 w0 = 1-w0;
101 102400 w1 = 1-w1;
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.
102400 PetscCall(MatCopy(C,ctx->W[w0],SAME_NONZERO_PATTERN));
103
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
102400 if (i == degree-1) {
104 512 s1 = s1/2;
105 512 s2 = s2/2;
106 }
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.
102400 PetscCall(MatScale(C,s2));
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.
102400 PetscCall(MatAXPY(C,-1,ctx->W[w1],SAME_NONZERO_PATTERN));
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.
102400 PetscCall(MatAXPY(C,cheby->damping_coeffs[degree-i-1]*cheby->coeffs[degree-i-1],B,SAME_NONZERO_PATTERN));
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.
102400 PetscCall(MatMatMult(ctx->T,ctx->W[w0],MAT_REUSE_MATRIX,PETSC_CURRENT,&ctx->W[w1]));
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.
102400 PetscCall(MatAXPY(C,s1,ctx->W[w1],SAME_NONZERO_PATTERN));
112 }
113
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.
73 PetscFunctionReturn(PETSC_SUCCESS);
114 }
115
116 98 static PetscErrorCode Chebyshev_Compute_Coefficients(ST st)
117 {
118 98 ST_FILTER *ctx = (ST_FILTER*)st->data;
119 98 CHEBYSHEV_CTX cheby = (CHEBYSHEV_CTX)ctx->data;
120 98 PetscInt i,degree;
121 98 PetscReal s_inta,s_intb,acosa,acosb,t1,t2,s;
122
123
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
98 PetscFunctionBegin;
124 98 degree = ctx->polyDegree;
125
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.
98 PetscCall(PetscFree(cheby->coeffs));
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.
98 PetscCall(PetscMalloc1(degree+1,&cheby->coeffs));
127
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
98 s_inta = PetscMax((2.0*ctx->inta - (ctx->right + ctx->left))/(ctx->right - ctx->left),-1.0);
128
1/2
✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
98 s_intb = PetscMin((2.0*ctx->intb - (ctx->right + ctx->left))/(ctx->right - ctx->left),1.0);
129 98 acosa = PetscAcosReal(s_inta);
130 98 acosb = PetscAcosReal(s_intb);
131 98 t1 = (acosa - acosb)/PETSC_PI;
132 98 t2 = t1;
133 98 cheby->coeffs[0] = t1;
134
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
13198 for (i=0;i<degree;i++) {
135 13100 cheby->coeffs[i+1] = 2.0*(PetscSinReal((i+1)*acosa) - PetscSinReal((i+1)*acosb)) / ((i+1)*PETSC_PI);
136 13100 t1 = t1 + cheby->coeffs[i+1]*PetscCosReal((i+1)*acosa);
137 13100 t2 = t2 + cheby->coeffs[i+1]*PetscCosReal((i+1)*acosb);
138 }
139
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
98 if (s_inta == -1.0) s = 0.5/PetscAbs(t1);
140
1/8
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
98 else if (s_intb == 1.0) s = 0.5/PetscMin(PetscAbs(t1),PetscAbs(t2));
141
4/6
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✓ Branch 4 taken 10 times.
✓ Branch 5 taken 10 times.
168 else s = 0.5/PetscMin(PetscAbs(t1),PetscAbs(t2));
142
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
13296 for (i=0;i<degree+1;i++) cheby->coeffs[i] = s*cheby->coeffs[i];
143
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.
18 PetscFunctionReturn(PETSC_SUCCESS);
144 }
145
146 98 static PetscErrorCode Chebyshev_Compute_Damping_Coefficients(ST st)
147 {
148 98 ST_FILTER *ctx = (ST_FILTER*)st->data;
149 98 CHEBYSHEV_CTX cheby = (CHEBYSHEV_CTX)ctx->data;
150 98 PetscInt i,degree;
151 98 PetscReal invdeg;
152
153
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
98 PetscFunctionBegin;
154 98 degree = ctx->polyDegree;
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.
98 PetscCall(PetscFree(cheby->damping_coeffs));
156
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.
98 PetscCall(PetscMalloc1(degree+1,&cheby->damping_coeffs));
157 98 cheby->damping_coeffs[0] = 1;
158
4/5
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
98 switch (ctx->damping) {
159 case ST_FILTER_DAMPING_NONE:
160
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
8148 for (i=1;i<=degree;i++) cheby->damping_coeffs[i] = 1;
161 break;
162 case ST_FILTER_DAMPING_LANCZOS:
163
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
2020 for (i=1;i<=degree;i++) cheby->damping_coeffs[i] = PetscPowRealInt(PetscSinReal(i*PETSC_PI/(degree+1))/(i*PETSC_PI/(degree+1)),DAMPING_LANCZOS_POW);
164 break;
165 10 case ST_FILTER_DAMPING_JACKSON:
166 10 invdeg = 1/((PetscReal)degree+2);
167
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
1010 for (i=1;i<=degree;i++) cheby->damping_coeffs[i] = invdeg*((degree+2-i)*PetscCosReal(PETSC_PI*i*invdeg)+PetscSinReal(PETSC_PI*i*invdeg)/PetscTanReal(PETSC_PI*invdeg));
168 break;
169 case ST_FILTER_DAMPING_FEJER:
170
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
2020 for (i=1;i<=degree;i++) cheby->damping_coeffs[i] = 1-i/(PetscReal)(degree+1);
171 break;
172 default:
173 SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_WRONG,"Invalid Chebyshev filter damping type");
174 }
175
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.
18 PetscFunctionReturn(PETSC_SUCCESS);
176 }
177
178 /*
179 Computes the coefficients and damping coefficients for the Chebyshev filter
180
181 Creates the shifted (and scaled) matrix and the filter P(z).
182 G is a shell matrix whose MatMult() and MatMatMult() applies the filter.
183 */
184 98 static PetscErrorCode STComputeOperator_Filter_Chebyshev(ST st,Mat *G)
185 {
186 98 ST_FILTER *ctx = (ST_FILTER*)st->data;
187 98 PetscInt n,m,N,M;
188
189
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
98 PetscFunctionBegin;
190
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.
98 PetscCall(STSetWorkVecs(st,2));
191
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.
98 PetscCall(Chebyshev_Compute_Coefficients(st));
192
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.
98 PetscCall(Chebyshev_Compute_Damping_Coefficients(st));
193 /* create shell matrix*/
194
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.
98 PetscCall(MatDestroy(&ctx->T));
195
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.
98 PetscCall(MatDuplicate(st->A[0],MAT_COPY_VALUES,&ctx->T));
196
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
98 if (!*G) {
197
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.
78 PetscCall(MatGetSize(ctx->T,&N,&M));
198
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.
78 PetscCall(MatGetLocalSize(ctx->T,&n,&m));
199
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.
78 PetscCall(MatCreateShell(PetscObjectComm((PetscObject)st),n,m,N,M,st,G));
200
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.
78 PetscCall(MatShellSetOperation(*G,MATOP_MULT,(PetscErrorCodeFn*)MatMult_Chebyshev));
201
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.
78 PetscCall(MatShellSetMatProductOperation(*G,MATPRODUCT_AB,NULL,MatMatMult_Chebyshev,NULL,MATDENSE,MATDENSE));
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.
78 PetscCall(MatShellSetMatProductOperation(*G,MATPRODUCT_AB,NULL,MatMatMult_Chebyshev,NULL,MATDENSECUDA,MATDENSECUDA));
203
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.
78 PetscCall(MatShellSetMatProductOperation(*G,MATPRODUCT_AB,NULL,MatMatMult_Chebyshev,NULL,MATDENSEHIP,MATDENSEHIP));
204 }
205
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.
18 PetscFunctionReturn(PETSC_SUCCESS);
206 }
207
208 334 static PetscErrorCode Chebyshev_Clenshaw_Damping_Real(ST st,PetscReal x,PetscReal *px)
209 {
210 334 ST_FILTER *ctx = (ST_FILTER*)st->data;
211 334 CHEBYSHEV_CTX cheby = (CHEBYSHEV_CTX)ctx->data;
212 334 PetscReal t1=0,t2=0;
213 334 PetscInt i,degree;
214
215
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
334 PetscFunctionBegin;
216 334 degree = ctx->polyDegree;
217
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
334 PetscCheck(cheby->coeffs,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_NULL,"Chebyshev coefficients are not computed");
218
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
334 PetscCheck(cheby->damping_coeffs,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_NULL,"Chebyshev damping coefficients are not computed");
219 334 x = 2*x;
220
2/2
✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
30234 for (i=degree;i>1;i=i-2) {
221 29900 t2 = cheby->damping_coeffs[i]*cheby->coeffs[i] + x*t1 - t2;
222 29900 t1 = cheby->damping_coeffs[i-1]*cheby->coeffs[i-1] + x*t2 - t1;
223 }
224
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
334 if (degree % 2) {
225 t2 = cheby->damping_coeffs[1]*cheby->coeffs[1] + x*t1 - t2;
226 t1 = t2;
227 }
228 334 *px = cheby->damping_coeffs[0]*cheby->coeffs[0] + .5*x*t1 - t2;
229
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.
334 PetscFunctionReturn(PETSC_SUCCESS);
230 }
231
232 /*
233 Computes the threshold for the Chebyshev filter
234
235 The threshold is P(scale(ctx->inta)), where ctx->inta is the lower bound of the interval.
236 It should be the same value if the upper bound was used instead.
237 */
238 334 static PetscErrorCode STFilterGetThreshold_Filter_Chebyshev(ST st,PetscReal *gamma)
239 {
240 334 ST_FILTER *ctx = (ST_FILTER*)st->data;
241 334 PetscReal c,e;
242
243
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
334 PetscFunctionBegin;
244 /* scale ctx->inta */
245 334 c = (ctx->right + ctx->left) / 2.0;
246 334 e = (ctx->right - ctx->left) / 2.0;
247
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.
334 PetscCall(Chebyshev_Clenshaw_Damping_Real(st,(ctx->inta-c)/e,gamma));
248
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.
52 PetscFunctionReturn(PETSC_SUCCESS);
249 }
250
251 78 static PetscErrorCode STDestroy_Filter_Chebyshev(ST st)
252 {
253 78 ST_FILTER *ctx = (ST_FILTER*)st->data;
254 78 CHEBYSHEV_CTX cheby = (CHEBYSHEV_CTX)ctx->data;
255
256
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
78 PetscFunctionBegin;
257
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.
78 PetscCall(PetscFree(cheby->coeffs));
258
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.
78 PetscCall(PetscFree(cheby->damping_coeffs));
259
6/8
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
78 PetscCall(PetscFree(cheby));
260
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.
14 PetscFunctionReturn(PETSC_SUCCESS);
261 }
262
263 78 PetscErrorCode STCreate_Filter_Chebyshev(ST st)
264 {
265 78 ST_FILTER *ctx = (ST_FILTER*)st->data;
266 78 CHEBYSHEV_CTX cheby;
267
268
1/2
✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
78 PetscFunctionBegin;
269
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.
78 PetscCall(PetscNew(&cheby));
270 78 ctx->data = (void*)cheby;
271
272 78 ctx->computeoperator = STComputeOperator_Filter_Chebyshev;
273 78 ctx->getthreshold = STFilterGetThreshold_Filter_Chebyshev;
274 78 ctx->destroy = STDestroy_Filter_Chebyshev;
275
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.
78 PetscFunctionReturn(PETSC_SUCCESS);
276 }
277