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[] = "Shell spectral transformations with a non-injective mapping.\n\n" | ||
12 | "Implements spectrum folding for the 2-D Laplacian, as in ex24.c.\n\n" | ||
13 | "The command line options are:\n" | ||
14 | " -n <n>, where <n> = number of grid subdivisions in x dimension.\n" | ||
15 | " -m <m>, where <m> = number of grid subdivisions in y dimension.\n"; | ||
16 | |||
17 | #include <slepceps.h> | ||
18 | |||
19 | /* Context for spectrum folding spectral transformation */ | ||
20 | typedef struct { | ||
21 | Mat A; | ||
22 | Vec w; | ||
23 | PetscScalar target; | ||
24 | } FoldShellST; | ||
25 | |||
26 | /* Routines for shell spectral transformation */ | ||
27 | PetscErrorCode STCreate_Fold(Mat,PetscScalar,FoldShellST**); | ||
28 | PetscErrorCode STApply_Fold(ST,Vec,Vec); | ||
29 | PetscErrorCode STDestroy_Fold(FoldShellST*); | ||
30 | |||
31 | 10 | int main (int argc,char **argv) | |
32 | { | ||
33 | 10 | Mat A; /* operator matrix */ | |
34 | 10 | EPS eps; /* eigenproblem solver context */ | |
35 | 10 | ST st; /* spectral transformation context */ | |
36 | 10 | FoldShellST *fold; /* user-defined spectral transform context */ | |
37 | 10 | EPSType type; | |
38 | 10 | PetscInt N,n=10,m,i,j,II,Istart,Iend,nev; | |
39 | 10 | PetscBool isShell,terse,flag; | |
40 | 10 | PetscScalar target=1.1; | |
41 | |||
42 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBeginUser; |
43 |
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(SlepcInitialize(&argc,&argv,NULL,help)); |
44 | |||
45 |
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(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); |
46 |
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(PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag)); |
47 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
10 | if (!flag) m = n; |
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.
|
10 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-target",&target,NULL)); |
49 | 10 | N = n*m; | |
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.
|
10 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nSpectrum Folding via shell ST, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid) target=%3.2f\n\n",N,n,m,(double)PetscRealPart(target))); |
51 | |||
52 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
53 | Compute the 5-point stencil Laplacian | ||
54 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
55 | |||
56 |
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(MatCreate(PETSC_COMM_WORLD,&A)); |
57 |
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(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
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.
|
10 | PetscCall(MatSetFromOptions(A)); |
59 | |||
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.
|
10 | PetscCall(MatGetOwnershipRange(A,&Istart,&Iend)); |
61 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1110 | for (II=Istart;II<Iend;II++) { |
62 | 1100 | i = II/n; j = II-i*n; | |
63 |
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.
|
1100 | if (i>0) PetscCall(MatSetValue(A,II,II-n,-1.0,INSERT_VALUES)); |
64 |
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.
|
1100 | if (i<m-1) PetscCall(MatSetValue(A,II,II+n,-1.0,INSERT_VALUES)); |
65 |
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.
|
1100 | if (j>0) PetscCall(MatSetValue(A,II,II-1,-1.0,INSERT_VALUES)); |
66 |
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.
|
1100 | if (j<n-1) PetscCall(MatSetValue(A,II,II+1,-1.0,INSERT_VALUES)); |
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.
|
1100 | PetscCall(MatSetValue(A,II,II,4.0,INSERT_VALUES)); |
68 | } | ||
69 |
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(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); |
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.
|
10 | PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); |
71 | |||
72 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
73 | Create the eigensolver and set various options | ||
74 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
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.
|
10 | PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps)); |
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.
|
10 | PetscCall(EPSSetOperators(eps,A,NULL)); |
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.
|
10 | PetscCall(EPSSetProblemType(eps,EPS_HEP)); |
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.
|
10 | PetscCall(EPSSetTarget(eps,target)); |
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.
|
10 | PetscCall(EPSGetST(eps,&st)); |
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.
|
10 | PetscCall(STSetType(st,STSHELL)); |
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.
|
10 | PetscCall(EPSSetFromOptions(eps)); |
83 | |||
84 | /* | ||
85 | Initialize shell spectral transformation | ||
86 | */ | ||
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.
|
10 | PetscCall(PetscObjectTypeCompare((PetscObject)st,STSHELL,&isShell)); |
88 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
10 | if (isShell) { |
89 | /* Change sorting criterion since this shell ST computes eigenvalues | ||
90 | of the transformed operator closest to 0 */ | ||
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.
|
10 | PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL)); |
92 | |||
93 | /* Create the context for the user-defined spectral transform */ | ||
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.
|
10 | PetscCall(STCreate_Fold(A,target,&fold)); |
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.
|
10 | PetscCall(STShellSetContext(st,fold)); |
96 | |||
97 | /* Set callback function for applying the operator (in this case we do not | ||
98 | provide a back-transformation callback since the mapping is not one-to-one) */ | ||
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.
|
10 | PetscCall(STShellSetApply(st,STApply_Fold)); |
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.
|
10 | PetscCall(PetscObjectSetName((PetscObject)st,"STFOLD")); |
101 | } | ||
102 | |||
103 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
104 | Solve the eigensystem | ||
105 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
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.
|
10 | PetscCall(EPSSolve(eps)); |
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.
|
10 | PetscCall(EPSGetType(eps,&type)); |
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.
|
10 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type)); |
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.
|
10 | PetscCall(EPSGetDimensions(eps,&nev,NULL,NULL)); |
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.
|
10 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of requested eigenvalues: %" PetscInt_FMT "\n",nev)); |
112 | |||
113 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
114 | Display solution and clean up | ||
115 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
116 | |||
117 | /* show detailed info unless -terse option is given by user */ | ||
118 |
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(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
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.
|
10 | if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL)); |
120 | else { | ||
121 | ✗ | PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL)); | |
122 | ✗ | PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD)); | |
123 | ✗ | PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD)); | |
124 | ✗ | PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); | |
125 | } | ||
126 |
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.
|
10 | if (isShell) PetscCall(STDestroy_Fold(fold)); |
127 |
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(EPSDestroy(&eps)); |
128 |
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(MatDestroy(&A)); |
129 |
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.
|
10 | PetscCall(SlepcFinalize()); |
130 | return 0; | ||
131 | } | ||
132 | |||
133 | /* | ||
134 | STCreate_Fold - Creates the spectrum folding ST context. | ||
135 | |||
136 | Input Parameter: | ||
137 | + A - problem matrix | ||
138 | - target - target value | ||
139 | |||
140 | Output Parameter: | ||
141 | . fold - user-defined spectral transformation context | ||
142 | */ | ||
143 | 10 | PetscErrorCode STCreate_Fold(Mat A,PetscScalar target,FoldShellST **fold) | |
144 | { | ||
145 | 10 | FoldShellST *newctx; | |
146 | |||
147 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBeginUser; |
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.
|
10 | PetscCall(PetscNew(&newctx)); |
149 | 10 | newctx->A = A; | |
150 |
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(PetscObjectReference((PetscObject)A)); |
151 | 10 | newctx->target = target; | |
152 |
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(MatCreateVecs(A,&newctx->w,NULL)); |
153 | 10 | *fold = newctx; | |
154 |
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.
|
10 | PetscFunctionReturn(PETSC_SUCCESS); |
155 | } | ||
156 | |||
157 | /* | ||
158 | STApply_Fold - Applies the operator (A-target*I)^2 to a given vector. | ||
159 | |||
160 | Input Parameters: | ||
161 | + st - spectral transformation context | ||
162 | - x - input vector | ||
163 | |||
164 | Output Parameter: | ||
165 | . y - output vector | ||
166 | */ | ||
167 | 5525 | PetscErrorCode STApply_Fold(ST st,Vec x,Vec y) | |
168 | { | ||
169 | 5525 | FoldShellST *fold; | |
170 | 5525 | PetscScalar sigma; | |
171 | |||
172 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
5525 | PetscFunctionBeginUser; |
173 |
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.
|
5525 | PetscCall(STShellGetContext(st,&fold)); |
174 | 5525 | sigma = -fold->target; | |
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.
|
5525 | PetscCall(MatMult(fold->A,x,fold->w)); |
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.
|
5525 | PetscCall(VecAXPY(fold->w,sigma,x)); |
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.
|
5525 | PetscCall(MatMult(fold->A,fold->w,y)); |
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.
|
5525 | PetscCall(VecAXPY(y,sigma,fold->w)); |
179 |
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.
|
1105 | PetscFunctionReturn(PETSC_SUCCESS); |
180 | } | ||
181 | |||
182 | /* | ||
183 | STDestroy_Fold - This routine destroys the shell ST context. | ||
184 | |||
185 | Input Parameter: | ||
186 | . fold - user-defined spectral transformation context | ||
187 | */ | ||
188 | 10 | PetscErrorCode STDestroy_Fold(FoldShellST *fold) | |
189 | { | ||
190 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBeginUser; |
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.
|
10 | PetscCall(MatDestroy(&fold->A)); |
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.
|
10 | PetscCall(VecDestroy(&fold->w)); |
193 |
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.
|
10 | PetscCall(PetscFree(fold)); |
194 |
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.
|
2 | PetscFunctionReturn(PETSC_SUCCESS); |
195 | } | ||
196 | |||
197 | /*TEST | ||
198 | |||
199 | test: | ||
200 | args: -m 11 -eps_nev 4 -terse | ||
201 | suffix: 1 | ||
202 | requires: !single | ||
203 | |||
204 | TEST*/ | ||
205 |