| 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[] = "Spectrum folding for a standard symmetric eigenproblem.\n\n" | ||
| 12 | "The problem matrix is the 2-D Laplacian.\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 | /* | ||
| 20 | User context for spectrum folding | ||
| 21 | */ | ||
| 22 | typedef struct { | ||
| 23 | Mat A; | ||
| 24 | Vec w; | ||
| 25 | PetscReal target; | ||
| 26 | } CTX_FOLD; | ||
| 27 | |||
| 28 | /* | ||
| 29 | Auxiliary routines | ||
| 30 | */ | ||
| 31 | PetscErrorCode MatMult_Fold(Mat,Vec,Vec); | ||
| 32 | PetscErrorCode RayleighQuotient(Mat,Vec,PetscScalar*); | ||
| 33 | PetscErrorCode ComputeResidualNorm(Mat,PetscScalar,Vec,PetscReal*); | ||
| 34 | |||
| 35 | 30 | int main(int argc,char **argv) | |
| 36 | { | ||
| 37 | 30 | Mat A,M,P; /* problem matrix, shell matrix and preconditioner */ | |
| 38 | 30 | Vec x; /* eigenvector */ | |
| 39 | 30 | EPS eps; /* eigenproblem solver context */ | |
| 40 | 30 | ST st; /* spectral transformation */ | |
| 41 | 30 | KSP ksp; | |
| 42 | 30 | PC pc; | |
| 43 | 30 | EPSType type; | |
| 44 | 30 | CTX_FOLD *ctx; | |
| 45 | 30 | PetscInt nconv,N,n=10,m,nloc,mloc,Istart,Iend,II,i,j; | |
| 46 | 30 | PetscReal *error,*evals,target=0.0,tol; | |
| 47 | 30 | PetscScalar lambda; | |
| 48 | 30 | PetscBool flag,terse,errok,hasmat; | |
| 49 | |||
| 50 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
30 | PetscFunctionBeginUser; |
| 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.
|
30 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
| 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.
|
30 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); |
| 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.
|
30 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag)); |
| 55 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
30 | if (!flag) m=n; |
| 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.
|
30 | PetscCall(PetscOptionsGetReal(NULL,NULL,"-target",&target,NULL)); |
| 57 | 30 | N = n*m; | |
| 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.
|
30 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nSpectrum Folding, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid) target=%f\n\n",N,n,m,(double)target)); |
| 59 | |||
| 60 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 61 | Compute the 5-point stencil Laplacian | ||
| 62 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 63 | |||
| 64 |
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(MatCreate(PETSC_COMM_WORLD,&A)); |
| 65 |
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(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
| 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.
|
30 | PetscCall(MatSetFromOptions(A)); |
| 67 | |||
| 68 |
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(MatGetOwnershipRange(A,&Istart,&Iend)); |
| 69 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
6780 | for (II=Istart;II<Iend;II++) { |
| 70 | 6750 | i = II/n; j = II-i*n; | |
| 71 |
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.
|
6750 | if (i>0) PetscCall(MatSetValue(A,II,II-n,-1.0,INSERT_VALUES)); |
| 72 |
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.
|
6750 | if (i<m-1) PetscCall(MatSetValue(A,II,II+n,-1.0,INSERT_VALUES)); |
| 73 |
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.
|
6750 | if (j>0) PetscCall(MatSetValue(A,II,II-1,-1.0,INSERT_VALUES)); |
| 74 |
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.
|
6750 | if (j<n-1) PetscCall(MatSetValue(A,II,II+1,-1.0,INSERT_VALUES)); |
| 75 |
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.
|
6750 | PetscCall(MatSetValue(A,II,II,4.0,INSERT_VALUES)); |
| 76 | } | ||
| 77 | |||
| 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.
|
30 | PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); |
| 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.
|
30 | PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); |
| 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.
|
30 | PetscCall(MatCreateVecs(A,&x,NULL)); |
| 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.
|
30 | PetscCall(MatGetLocalSize(A,&nloc,&mloc)); |
| 82 | |||
| 83 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 84 | Create shell matrix to perform spectrum folding | ||
| 85 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 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.
|
30 | PetscCall(PetscNew(&ctx)); |
| 87 | 30 | ctx->A = A; | |
| 88 | 30 | ctx->target = target; | |
| 89 |
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(VecDuplicate(x,&ctx->w)); |
| 90 | |||
| 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.
|
30 | PetscCall(MatCreateShell(PETSC_COMM_WORLD,nloc,mloc,N,N,ctx,&M)); |
| 92 |
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(MatShellSetOperation(M,MATOP_MULT,(PetscErrorCodeFn*)MatMult_Fold)); |
| 93 | |||
| 94 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 95 | Create the eigensolver and set various options | ||
| 96 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 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.
|
30 | PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps)); |
| 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.
|
30 | PetscCall(EPSSetOperators(eps,M,NULL)); |
| 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.
|
30 | PetscCall(EPSSetProblemType(eps,EPS_HEP)); |
| 101 |
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_SMALLEST_REAL)); |
| 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.
|
30 | PetscCall(EPSSetFromOptions(eps)); |
| 103 | |||
| 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.
|
30 | PetscCall(PetscObjectTypeCompareAny((PetscObject)eps,&flag,EPSGD,EPSJD,EPSBLOPEX,EPSLOBPCG,EPSRQCG,"")); |
| 105 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
30 | if (flag) { |
| 106 | /* | ||
| 107 | Build preconditioner specific for this application (diagonal of A^2) | ||
| 108 | */ | ||
| 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.
|
20 | PetscCall(MatGetRowSum(A,x)); |
| 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.
|
20 | PetscCall(VecScale(x,-1.0)); |
| 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.
|
20 | PetscCall(VecShift(x,20.0)); |
| 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.
|
20 | PetscCall(MatCreate(PETSC_COMM_WORLD,&P)); |
| 113 |
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.
|
20 | PetscCall(MatSetSizes(P,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
| 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.
|
20 | PetscCall(MatSetFromOptions(P)); |
| 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.
|
20 | PetscCall(MatDiagonalSet(P,x,INSERT_VALUES)); |
| 116 | /* | ||
| 117 | Set diagonal preconditioner | ||
| 118 | */ | ||
| 119 |
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.
|
20 | PetscCall(EPSGetST(eps,&st)); |
| 120 |
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.
|
20 | PetscCall(STSetType(st,STPRECOND)); |
| 121 |
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.
|
20 | PetscCall(STSetPreconditionerMat(st,P)); |
| 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.
|
20 | PetscCall(MatDestroy(&P)); |
| 123 |
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.
|
20 | PetscCall(STGetKSP(st,&ksp)); |
| 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.
|
20 | PetscCall(KSPGetPC(ksp,&pc)); |
| 125 |
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.
|
20 | PetscCall(PCSetType(pc,PCJACOBI)); |
| 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.
|
20 | PetscCall(STPrecondGetKSPHasMat(st,&hasmat)); |
| 127 |
7/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
30 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Preconditioned solver, hasmat=%s\n",hasmat?"true":"false")); |
| 128 | } | ||
| 129 | |||
| 130 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 131 | Solve the eigensystem | ||
| 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.
|
30 | PetscCall(EPSSolve(eps)); |
| 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.
|
30 | PetscCall(EPSGetType(eps,&type)); |
| 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.
|
30 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type)); |
| 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.
|
30 | PetscCall(EPSGetTolerances(eps,&tol,NULL)); |
| 138 | |||
| 139 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 140 | Display solution and clean up | ||
| 141 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 142 | |||
| 143 |
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(EPSGetConverged(eps,&nconv)); |
| 144 |
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(PetscPrintf(PETSC_COMM_WORLD," Number of converged eigenpairs: %" PetscInt_FMT "\n\n",nconv)); |
| 145 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
30 | if (nconv>0) { |
| 146 |
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(PetscMalloc2(nconv,&evals,nconv,&error)); |
| 147 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
60 | for (i=0;i<nconv;i++) { |
| 148 | /* Get i-th eigenvector, compute eigenvalue approximation from | ||
| 149 | Rayleigh quotient and compute residual norm */ | ||
| 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.
|
30 | PetscCall(EPSGetEigenpair(eps,i,NULL,NULL,x,NULL)); |
| 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.
|
30 | PetscCall(RayleighQuotient(A,x,&lambda)); |
| 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.
|
30 | PetscCall(ComputeResidualNorm(A,lambda,x,&error[i])); |
| 153 | #if defined(PETSC_USE_COMPLEX) | ||
| 154 | 15 | evals[i] = PetscRealPart(lambda); | |
| 155 | #else | ||
| 156 | 15 | evals[i] = lambda; | |
| 157 | #endif | ||
| 158 | } | ||
| 159 |
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(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
| 160 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
30 | if (!terse) { |
| 161 | ✗ | PetscCall(PetscPrintf(PETSC_COMM_WORLD, | |
| 162 | " k ||Ax-kx||\n" | ||
| 163 | " ----------------- ------------------\n")); | ||
| 164 | ✗ | for (i=0;i<nconv;i++) PetscCall(PetscPrintf(PETSC_COMM_WORLD," %12f %12.2g\n",(double)evals[i],(double)error[i])); | |
| 165 | } else { | ||
| 166 | errok = PETSC_TRUE; | ||
| 167 |
4/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✓ Branch 4 taken 10 times.
✓ Branch 5 taken 10 times.
|
60 | for (i=0;i<nconv;i++) errok = (errok && error[i]<5.0*tol)? PETSC_TRUE: PETSC_FALSE; |
| 168 |
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.
|
30 | if (!errok) PetscCall(PetscPrintf(PETSC_COMM_WORLD," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nconv)); |
| 169 | else { | ||
| 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.
|
30 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," nconv=%" PetscInt_FMT " eigenvalues computed up to the required tolerance:",nconv)); |
| 171 |
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.
|
60 | for (i=0;i<nconv;i++) PetscCall(PetscPrintf(PETSC_COMM_WORLD," %.5f",(double)evals[i])); |
| 172 | } | ||
| 173 | } | ||
| 174 |
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(PetscPrintf(PETSC_COMM_WORLD,"\n")); |
| 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.
|
30 | PetscCall(PetscFree2(evals,error)); |
| 176 | } | ||
| 177 | |||
| 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.
|
30 | PetscCall(EPSDestroy(&eps)); |
| 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.
|
30 | PetscCall(MatDestroy(&A)); |
| 180 |
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(MatDestroy(&M)); |
| 181 |
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(VecDestroy(&ctx->w)); |
| 182 |
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(VecDestroy(&x)); |
| 183 |
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.
|
30 | PetscCall(PetscFree(ctx)); |
| 184 |
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.
|
30 | PetscCall(SlepcFinalize()); |
| 185 | return 0; | ||
| 186 | } | ||
| 187 | |||
| 188 | /* | ||
| 189 | Matrix-vector product subroutine for the spectrum folding. | ||
| 190 | y <-- (A-t*I)^2*x | ||
| 191 | */ | ||
| 192 | 13170 | PetscErrorCode MatMult_Fold(Mat M,Vec x,Vec y) | |
| 193 | { | ||
| 194 | 13170 | CTX_FOLD *ctx; | |
| 195 | 13170 | PetscScalar sigma; | |
| 196 | |||
| 197 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
13170 | PetscFunctionBeginUser; |
| 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.
|
13170 | PetscCall(MatShellGetContext(M,&ctx)); |
| 199 | 13170 | sigma = -ctx->target; | |
| 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.
|
13170 | PetscCall(MatMult(ctx->A,x,ctx->w)); |
| 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.
|
13170 | PetscCall(VecAXPY(ctx->w,sigma,x)); |
| 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.
|
13170 | PetscCall(MatMult(ctx->A,ctx->w,y)); |
| 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.
|
13170 | PetscCall(VecAXPY(y,sigma,ctx->w)); |
| 204 |
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.
|
2634 | PetscFunctionReturn(PETSC_SUCCESS); |
| 205 | } | ||
| 206 | |||
| 207 | /* | ||
| 208 | Computes the Rayleigh quotient of a vector x | ||
| 209 | r <-- x^T*A*x (assumes x has unit norm) | ||
| 210 | */ | ||
| 211 | 30 | PetscErrorCode RayleighQuotient(Mat A,Vec x,PetscScalar *r) | |
| 212 | { | ||
| 213 | 30 | Vec Ax; | |
| 214 | |||
| 215 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
30 | PetscFunctionBeginUser; |
| 216 |
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(VecDuplicate(x,&Ax)); |
| 217 |
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(MatMult(A,x,Ax)); |
| 218 |
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(VecDot(Ax,x,r)); |
| 219 |
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(VecDestroy(&Ax)); |
| 220 |
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.
|
6 | PetscFunctionReturn(PETSC_SUCCESS); |
| 221 | } | ||
| 222 | |||
| 223 | /* | ||
| 224 | Computes the residual norm of an approximate eigenvector x, |A*x-lambda*x| | ||
| 225 | */ | ||
| 226 | 30 | PetscErrorCode ComputeResidualNorm(Mat A,PetscScalar lambda,Vec x,PetscReal *r) | |
| 227 | { | ||
| 228 | 30 | Vec Ax; | |
| 229 | |||
| 230 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
30 | PetscFunctionBeginUser; |
| 231 |
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(VecDuplicate(x,&Ax)); |
| 232 |
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(MatMult(A,x,Ax)); |
| 233 |
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(VecAXPY(Ax,-lambda,x)); |
| 234 |
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(VecNorm(Ax,NORM_2,r)); |
| 235 |
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(VecDestroy(&Ax)); |
| 236 |
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.
|
6 | PetscFunctionReturn(PETSC_SUCCESS); |
| 237 | } | ||
| 238 | |||
| 239 | /*TEST | ||
| 240 | |||
| 241 | testset: | ||
| 242 | args: -n 15 -eps_nev 1 -eps_ncv 12 -eps_max_it 1000 -eps_tol 1e-5 -terse | ||
| 243 | filter: grep -v Solution | ||
| 244 | test: | ||
| 245 | suffix: 1 | ||
| 246 | test: | ||
| 247 | suffix: 1_lobpcg | ||
| 248 | args: -eps_type lobpcg | ||
| 249 | requires: !single | ||
| 250 | test: | ||
| 251 | suffix: 1_gd | ||
| 252 | args: -eps_type gd | ||
| 253 | requires: !single | ||
| 254 | |||
| 255 | TEST*/ | ||
| 256 |