| 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 | SLEPc quadratic eigensolver: "qarnoldi" | ||
| 12 | |||
| 13 | Method: Q-Arnoldi | ||
| 14 | |||
| 15 | Algorithm: | ||
| 16 | |||
| 17 | Quadratic Arnoldi with Krylov-Schur type restart. | ||
| 18 | |||
| 19 | References: | ||
| 20 | |||
| 21 | [1] K. Meerbergen, "The Quadratic Arnoldi method for the solution | ||
| 22 | of the quadratic eigenvalue problem", SIAM J. Matrix Anal. | ||
| 23 | Appl. 30(4):1462-1482, 2009. | ||
| 24 | */ | ||
| 25 | |||
| 26 | #include <slepc/private/pepimpl.h> /*I "slepcpep.h" I*/ | ||
| 27 | #include <petscblaslapack.h> | ||
| 28 | |||
| 29 | typedef struct { | ||
| 30 | PetscReal keep; /* restart parameter */ | ||
| 31 | PetscBool lock; /* locking/non-locking variant */ | ||
| 32 | } PEP_QARNOLDI; | ||
| 33 | |||
| 34 | 179 | static PetscErrorCode PEPSetUp_QArnoldi(PEP pep) | |
| 35 | { | ||
| 36 | 179 | PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data; | |
| 37 | 179 | PetscBool flg; | |
| 38 | |||
| 39 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
179 | PetscFunctionBegin; |
| 40 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
179 | PEPCheckQuadratic(pep); |
| 41 |
6/10✓ 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.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
|
179 | PEPCheckShiftSinvert(pep); |
| 42 |
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.
|
179 | PetscCall(PEPSetDimensions_Default(pep,pep->nev,&pep->ncv,&pep->mpd)); |
| 43 |
3/6✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
179 | PetscCheck(ctx->lock || pep->mpd>=pep->ncv,PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Should not use mpd parameter in non-locking variant"); |
| 44 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
179 | if (pep->max_it==PETSC_DETERMINE) pep->max_it = PetscMax(100,4*pep->n/pep->ncv); |
| 45 |
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.
|
179 | if (!pep->which) PetscCall(PEPSetWhichEigenpairs_Default(pep)); |
| 46 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
179 | PetscCheck(pep->which!=PEP_ALL,PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"This solver does not support computing all eigenvalues"); |
| 47 | |||
| 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.
|
179 | PetscCall(STGetTransform(pep->st,&flg)); |
| 49 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
179 | PetscCheck(flg,PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Solver requires the ST transformation flag set, see STSetTransform()"); |
| 50 | |||
| 51 | /* set default extraction */ | ||
| 52 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
179 | if (!pep->extract) pep->extract = PEP_EXTRACT_NONE; |
| 53 |
2/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 10 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
179 | PEPCheckUnsupported(pep,PEP_FEATURE_NONMONOMIAL | PEP_FEATURE_EXTRACT); |
| 54 | |||
| 55 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
179 | if (!ctx->keep) ctx->keep = 0.5; |
| 56 | |||
| 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.
|
179 | PetscCall(PEPAllocateSolution(pep,0)); |
| 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.
|
179 | PetscCall(PEPSetWorkVecs(pep,4)); |
| 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.
|
179 | PetscCall(DSSetType(pep->ds,DSNHEP)); |
| 61 |
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.
|
179 | PetscCall(DSSetExtraRow(pep->ds,PETSC_TRUE)); |
| 62 |
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.
|
179 | PetscCall(DSAllocate(pep->ds,pep->ncv+1)); |
| 63 |
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.
|
35 | PetscFunctionReturn(PETSC_SUCCESS); |
| 64 | } | ||
| 65 | |||
| 66 | 179 | static PetscErrorCode PEPExtractVectors_QArnoldi(PEP pep) | |
| 67 | { | ||
| 68 | 179 | PetscInt k=pep->nconv; | |
| 69 | 179 | Mat X,X0; | |
| 70 | |||
| 71 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
179 | PetscFunctionBegin; |
| 72 |
2/14✓ Branch 0 taken 8 times.
✓ Branch 1 taken 2 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.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
179 | if (pep->nconv==0) PetscFunctionReturn(PETSC_SUCCESS); |
| 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.
|
179 | PetscCall(DSVectors(pep->ds,DS_MAT_X,NULL,NULL)); |
| 74 | |||
| 75 | /* update vectors V = V*X */ | ||
| 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.
|
179 | PetscCall(DSGetMat(pep->ds,DS_MAT_X,&X)); |
| 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.
|
179 | PetscCall(MatDenseGetSubMatrix(X,0,k,0,k,&X0)); |
| 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.
|
179 | PetscCall(BVMultInPlace(pep->V,X0,0,k)); |
| 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.
|
179 | PetscCall(MatDenseRestoreSubMatrix(X,&X0)); |
| 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.
|
179 | PetscCall(DSRestoreMat(pep->ds,DS_MAT_X,&X)); |
| 81 |
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.
|
35 | PetscFunctionReturn(PETSC_SUCCESS); |
| 82 | } | ||
| 83 | |||
| 84 | /* | ||
| 85 | Compute a step of Classical Gram-Schmidt orthogonalization | ||
| 86 | */ | ||
| 87 | 29807 | static PetscErrorCode PEPQArnoldiCGS(PEP pep,PetscScalar *H,PetscBLASInt ldh,PetscScalar *h,PetscBLASInt j,BV V,Vec t,Vec v,Vec w,PetscReal *onorm,PetscReal *norm,PetscScalar *work) | |
| 88 | { | ||
| 89 | 29807 | PetscBLASInt ione = 1,j_1 = j+1; | |
| 90 | 29807 | PetscReal x,y; | |
| 91 | 29807 | PetscScalar dot,one = 1.0,zero = 0.0; | |
| 92 | |||
| 93 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
29807 | PetscFunctionBegin; |
| 94 | /* compute norm of v and w */ | ||
| 95 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
29807 | if (onorm) { |
| 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.
|
14897 | PetscCall(VecNorm(v,NORM_2,&x)); |
| 97 |
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.
|
14897 | PetscCall(VecNorm(w,NORM_2,&y)); |
| 98 | 14897 | *onorm = SlepcAbs(x,y); | |
| 99 | } | ||
| 100 | |||
| 101 | /* orthogonalize: compute h */ | ||
| 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.
|
29807 | PetscCall(BVDotVec(V,v,h)); |
| 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.
|
29807 | PetscCall(BVDotVec(V,w,work)); |
| 104 |
12/22✓ 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 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✓ Branch 19 taken 2 times.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
|
29807 | if (j>0) PetscCallBLAS("BLASgemv",BLASgemv_("C",&j_1,&j,&one,H,&ldh,work,&ione,&one,h,&ione)); |
| 105 |
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.
|
29807 | PetscCall(VecDot(w,t,&dot)); |
| 106 | 29807 | h[j] += dot; | |
| 107 | |||
| 108 | /* orthogonalize: update v and w */ | ||
| 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.
|
29807 | PetscCall(BVMultVec(V,-1.0,1.0,v,h)); |
| 110 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
29807 | if (j>0) { |
| 111 |
10/20✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✓ 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 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
|
29608 | PetscCallBLAS("BLASgemv",BLASgemv_("N",&j_1,&j,&one,H,&ldh,h,&ione,&zero,work,&ione)); |
| 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.
|
29608 | PetscCall(BVMultVec(V,-1.0,1.0,w,work)); |
| 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.
|
29807 | PetscCall(VecAXPY(w,-h[j],t)); |
| 115 | |||
| 116 | /* compute norm of v and w */ | ||
| 117 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
29807 | if (norm) { |
| 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.
|
28512 | PetscCall(VecNorm(v,NORM_2,&x)); |
| 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.
|
28512 | PetscCall(VecNorm(w,NORM_2,&y)); |
| 120 | 28512 | *norm = SlepcAbs(x,y); | |
| 121 | } | ||
| 122 |
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.
|
5913 | PetscFunctionReturn(PETSC_SUCCESS); |
| 123 | } | ||
| 124 | |||
| 125 | /* | ||
| 126 | Compute a run of Q-Arnoldi iterations | ||
| 127 | */ | ||
| 128 | 1270 | static PetscErrorCode PEPQArnoldi(PEP pep,Mat A,PetscInt k,PetscInt *M,Vec v,Vec w,PetscReal *beta,PetscBool *breakdown,PetscScalar *work) | |
| 129 | { | ||
| 130 | 1270 | PetscInt i,j,l,m = *M,ldh; | |
| 131 | 1270 | PetscBLASInt jj,ldhh; | |
| 132 | 1270 | Vec t = pep->work[2],u = pep->work[3]; | |
| 133 | 1270 | BVOrthogRefineType refinement; | |
| 134 | 1270 | PetscReal norm=0.0,onorm,eta; | |
| 135 | 1270 | PetscScalar *H,*c = work + m; | |
| 136 | |||
| 137 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
1270 | PetscFunctionBegin; |
| 138 | 1270 | *beta = 0.0; | |
| 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.
|
1270 | PetscCall(MatDenseGetArray(A,&H)); |
| 140 |
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.
|
1270 | PetscCall(MatDenseGetLDA(A,&ldh)); |
| 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.
|
1270 | PetscCall(BVGetOrthogonalization(pep->V,NULL,&refinement,&eta,NULL)); |
| 142 |
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.
|
1270 | PetscCall(BVInsertVec(pep->V,k,v)); |
| 143 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
16993 | for (j=k;j<m;j++) { |
| 144 | /* apply operator */ | ||
| 145 |
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.
|
15723 | PetscCall(VecCopy(w,t)); |
| 146 |
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.
|
15723 | if (pep->Dr) PetscCall(VecPointwiseMult(v,v,pep->Dr)); |
| 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.
|
15723 | PetscCall(STMatMult(pep->st,0,v,u)); |
| 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.
|
15723 | PetscCall(VecCopy(t,v)); |
| 149 |
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.
|
15723 | if (pep->Dr) PetscCall(VecPointwiseMult(t,t,pep->Dr)); |
| 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.
|
15723 | PetscCall(STMatMult(pep->st,1,t,w)); |
| 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.
|
15723 | PetscCall(VecAXPY(u,pep->sfactor,w)); |
| 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.
|
15723 | PetscCall(STMatSolve(pep->st,u,w)); |
| 153 |
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.
|
15723 | PetscCall(VecScale(w,-1.0/(pep->sfactor*pep->sfactor))); |
| 154 |
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.
|
15723 | if (pep->Dr) PetscCall(VecPointwiseDivide(w,w,pep->Dr)); |
| 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.
|
15723 | PetscCall(VecCopy(v,t)); |
| 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.
|
15723 | PetscCall(BVSetActiveColumns(pep->V,0,j+1)); |
| 157 | |||
| 158 | /* orthogonalize */ | ||
| 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.
|
15723 | PetscCall(PetscBLASIntCast(j,&jj)); |
| 160 |
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.
|
15723 | PetscCall(PetscBLASIntCast(ldh,&ldhh)); |
| 161 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
15723 | switch (refinement) { |
| 162 | 826 | case BV_ORTHOG_REFINE_NEVER: | |
| 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.
|
826 | PetscCall(PEPQArnoldiCGS(pep,H,ldhh,H+ldh*j,jj,pep->V,t,v,w,NULL,&norm,work)); |
| 164 | 826 | *breakdown = PETSC_FALSE; | |
| 165 | 826 | break; | |
| 166 | 1295 | case BV_ORTHOG_REFINE_ALWAYS: | |
| 167 |
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.
|
1295 | PetscCall(PEPQArnoldiCGS(pep,H,ldhh,H+ldh*j,jj,pep->V,t,v,w,NULL,NULL,work)); |
| 168 |
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.
|
1295 | PetscCall(PEPQArnoldiCGS(pep,H,ldhh,c,jj,pep->V,t,v,w,&onorm,&norm,work)); |
| 169 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
20390 | for (i=0;i<=j;i++) H[ldh*j+i] += c[i]; |
| 170 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
1295 | if (norm < eta * onorm) *breakdown = PETSC_TRUE; |
| 171 | 1295 | else *breakdown = PETSC_FALSE; | |
| 172 | break; | ||
| 173 | 13602 | case BV_ORTHOG_REFINE_IFNEEDED: | |
| 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.
|
13602 | PetscCall(PEPQArnoldiCGS(pep,H,ldhh,H+ldh*j,jj,pep->V,t,v,w,&onorm,&norm,work)); |
| 175 | /* ||q|| < eta ||h|| */ | ||
| 176 | l = 1; | ||
| 177 |
3/4✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
26391 | while (l<3 && norm < eta * onorm) { |
| 178 | 12789 | l++; | |
| 179 | 12789 | onorm = norm; | |
| 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.
|
12789 | PetscCall(PEPQArnoldiCGS(pep,H,ldhh,c,jj,pep->V,t,v,w,NULL,&norm,work)); |
| 181 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
231725 | for (i=0;i<=j;i++) H[ldh*j+i] += c[i]; |
| 182 | } | ||
| 183 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
13602 | if (norm < eta * onorm) *breakdown = PETSC_TRUE; |
| 184 | 13602 | else *breakdown = PETSC_FALSE; | |
| 185 | break; | ||
| 186 | } | ||
| 187 |
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.
|
15723 | PetscCall(VecScale(v,1.0/norm)); |
| 188 |
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.
|
15723 | PetscCall(VecScale(w,1.0/norm)); |
| 189 | |||
| 190 | 15723 | H[j+1+ldh*j] = norm; | |
| 191 |
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.
|
15723 | if (j<m-1) PetscCall(BVInsertVec(pep->V,j+1,v)); |
| 192 | } | ||
| 193 | 1270 | *beta = norm; | |
| 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.
|
1270 | PetscCall(MatDenseRestoreArray(A,&H)); |
| 195 |
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.
|
250 | PetscFunctionReturn(PETSC_SUCCESS); |
| 196 | } | ||
| 197 | |||
| 198 | 179 | static PetscErrorCode PEPSolve_QArnoldi(PEP pep) | |
| 199 | { | ||
| 200 | 179 | PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data; | |
| 201 | 179 | PetscInt j,k,l,lwork,nv,nconv; | |
| 202 | 179 | Vec v=pep->work[0],w=pep->work[1]; | |
| 203 | 179 | Mat Q,S; | |
| 204 | 179 | PetscScalar *work; | |
| 205 | 179 | PetscReal beta,norm,x,y; | |
| 206 | 179 | PetscBool breakdown=PETSC_FALSE,sinv; | |
| 207 | |||
| 208 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
179 | PetscFunctionBegin; |
| 209 | 179 | lwork = 7*pep->ncv; | |
| 210 |
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.
|
179 | PetscCall(PetscMalloc1(lwork,&work)); |
| 211 |
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.
|
179 | PetscCall(PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv)); |
| 212 |
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.
|
179 | PetscCall(RGPushScale(pep->rg,sinv?pep->sfactor:1.0/pep->sfactor)); |
| 213 |
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.
|
179 | PetscCall(STScaleShift(pep->st,sinv?pep->sfactor:1.0/pep->sfactor)); |
| 214 | |||
| 215 | /* Get the starting Arnoldi vector */ | ||
| 216 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
537 | for (j=0;j<2;j++) { |
| 217 |
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.
|
358 | if (j>=pep->nini) PetscCall(BVSetRandomColumn(pep->V,j)); |
| 218 | } | ||
| 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.
|
179 | PetscCall(BVCopyVec(pep->V,0,v)); |
| 220 |
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.
|
179 | PetscCall(BVCopyVec(pep->V,1,w)); |
| 221 |
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.
|
179 | PetscCall(VecNorm(v,NORM_2,&x)); |
| 222 |
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.
|
179 | PetscCall(VecNorm(w,NORM_2,&y)); |
| 223 | 179 | norm = SlepcAbs(x,y); | |
| 224 |
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.
|
179 | PetscCall(VecScale(v,1.0/norm)); |
| 225 |
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.
|
179 | PetscCall(VecScale(w,1.0/norm)); |
| 226 | |||
| 227 | /* clean projected matrix (including the extra-arrow) */ | ||
| 228 |
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.
|
179 | PetscCall(DSSetDimensions(pep->ds,PETSC_DETERMINE,PETSC_DETERMINE,PETSC_DETERMINE)); |
| 229 |
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.
|
179 | PetscCall(DSGetMat(pep->ds,DS_MAT_A,&S)); |
| 230 |
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.
|
179 | PetscCall(MatZeroEntries(S)); |
| 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.
|
179 | PetscCall(DSRestoreMat(pep->ds,DS_MAT_A,&S)); |
| 232 | |||
| 233 | /* Restart loop */ | ||
| 234 | 179 | l = 0; | |
| 235 | 179 | while (pep->reason == PEP_CONVERGED_ITERATING) { | |
| 236 | 1270 | pep->its++; | |
| 237 | |||
| 238 | /* Compute an nv-step Arnoldi factorization */ | ||
| 239 | 1270 | nv = PetscMin(pep->nconv+pep->mpd,pep->ncv); | |
| 240 |
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.
|
1270 | PetscCall(DSGetMat(pep->ds,DS_MAT_A,&S)); |
| 241 |
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.
|
1270 | PetscCall(PEPQArnoldi(pep,S,pep->nconv+l,&nv,v,w,&beta,&breakdown,work)); |
| 242 |
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.
|
1270 | PetscCall(DSRestoreMat(pep->ds,DS_MAT_A,&S)); |
| 243 |
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.
|
1270 | PetscCall(DSSetDimensions(pep->ds,nv,pep->nconv,pep->nconv+l)); |
| 244 |
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.
|
1270 | PetscCall(DSSetState(pep->ds,l?DS_STATE_RAW:DS_STATE_INTERMEDIATE)); |
| 245 |
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.
|
1270 | PetscCall(BVSetActiveColumns(pep->V,pep->nconv,nv)); |
| 246 | |||
| 247 | /* Solve projected problem */ | ||
| 248 |
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.
|
1270 | PetscCall(DSSolve(pep->ds,pep->eigr,pep->eigi)); |
| 249 |
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.
|
1270 | PetscCall(DSSort(pep->ds,pep->eigr,pep->eigi,NULL,NULL,NULL)); |
| 250 |
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.
|
1270 | PetscCall(DSUpdateExtraRow(pep->ds)); |
| 251 |
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.
|
1270 | PetscCall(DSSynchronize(pep->ds,pep->eigr,pep->eigi)); |
| 252 | |||
| 253 | /* Check convergence */ | ||
| 254 |
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.
|
1270 | PetscCall(PEPKrylovConvergence(pep,PETSC_FALSE,pep->nconv,nv-pep->nconv,beta,&k)); |
| 255 |
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.
|
1270 | PetscCall((*pep->stopping)(pep,pep->its,pep->max_it,k,pep->nev,&pep->reason,pep->stoppingctx)); |
| 256 | 1270 | nconv = k; | |
| 257 | |||
| 258 | /* Update l */ | ||
| 259 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
1270 | if (pep->reason != PEP_CONVERGED_ITERATING || breakdown) l = 0; |
| 260 | else { | ||
| 261 | 1091 | l = PetscMax(1,(PetscInt)((nv-k)*ctx->keep)); | |
| 262 |
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.
|
1091 | PetscCall(DSGetTruncateSize(pep->ds,k,nv,&l)); |
| 263 | } | ||
| 264 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
1270 | if (!ctx->lock && l>0) { l += k; k = 0; } /* non-locking variant: reset no. of converged pairs */ |
| 265 |
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.
|
1270 | if (l) PetscCall(PetscInfo(pep,"Preparing to restart keeping l=%" PetscInt_FMT " vectors\n",l)); |
| 266 | |||
| 267 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1270 | if (pep->reason == PEP_CONVERGED_ITERATING) { |
| 268 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
1091 | if (PetscUnlikely(breakdown)) { |
| 269 | /* Stop if breakdown */ | ||
| 270 | ✗ | PetscCall(PetscInfo(pep,"Breakdown Quadratic Arnoldi method (it=%" PetscInt_FMT " norm=%g)\n",pep->its,(double)beta)); | |
| 271 | ✗ | pep->reason = PEP_DIVERGED_BREAKDOWN; | |
| 272 | } else { | ||
| 273 | /* Prepare the Rayleigh quotient for restart */ | ||
| 274 |
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.
|
1091 | PetscCall(DSTruncate(pep->ds,k+l,PETSC_FALSE)); |
| 275 | } | ||
| 276 | } | ||
| 277 | /* Update the corresponding vectors V(:,idx) = V*Q(:,idx) */ | ||
| 278 |
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.
|
1270 | PetscCall(DSGetMat(pep->ds,DS_MAT_Q,&Q)); |
| 279 |
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.
|
1270 | PetscCall(BVMultInPlace(pep->V,Q,pep->nconv,k+l)); |
| 280 |
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.
|
1270 | PetscCall(DSRestoreMat(pep->ds,DS_MAT_Q,&Q)); |
| 281 | |||
| 282 | 1270 | pep->nconv = k; | |
| 283 |
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.
|
1449 | PetscCall(PEPMonitor(pep,pep->its,nconv,pep->eigr,pep->eigi,pep->errest,nv)); |
| 284 | } | ||
| 285 |
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.
|
179 | PetscCall(BVSetActiveColumns(pep->V,0,pep->nconv)); |
| 286 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1137 | for (j=0;j<pep->nconv;j++) { |
| 287 | 958 | pep->eigr[j] *= pep->sfactor; | |
| 288 | 958 | pep->eigi[j] *= pep->sfactor; | |
| 289 | } | ||
| 290 | |||
| 291 |
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.
|
179 | PetscCall(STScaleShift(pep->st,sinv?1.0/pep->sfactor:pep->sfactor)); |
| 292 |
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.
|
179 | PetscCall(RGPopScale(pep->rg)); |
| 293 | |||
| 294 |
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.
|
179 | PetscCall(DSTruncate(pep->ds,pep->nconv,PETSC_TRUE)); |
| 295 |
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.
|
179 | PetscCall(PetscFree(work)); |
| 296 |
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.
|
35 | PetscFunctionReturn(PETSC_SUCCESS); |
| 297 | } | ||
| 298 | |||
| 299 | 10 | static PetscErrorCode PEPQArnoldiSetRestart_QArnoldi(PEP pep,PetscReal keep) | |
| 300 | { | ||
| 301 | 10 | PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data; | |
| 302 | |||
| 303 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBegin; |
| 304 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
10 | if (keep==(PetscReal)PETSC_DEFAULT || keep==(PetscReal)PETSC_DECIDE) ctx->keep = 0.5; |
| 305 | else { | ||
| 306 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
10 | PetscCheck(keep>=0.1 && keep<=0.9,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"The keep argument must be in the range [0.1,0.9]"); |
| 307 | 10 | ctx->keep = keep; | |
| 308 | } | ||
| 309 |
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.
|
2 | PetscFunctionReturn(PETSC_SUCCESS); |
| 310 | } | ||
| 311 | |||
| 312 | /*@ | ||
| 313 | PEPQArnoldiSetRestart - Sets the restart parameter for the Q-Arnoldi | ||
| 314 | method, in particular the proportion of basis vectors that must be kept | ||
| 315 | after restart. | ||
| 316 | |||
| 317 | Logically Collective | ||
| 318 | |||
| 319 | Input Parameters: | ||
| 320 | + pep - the polynomial eigensolver context | ||
| 321 | - keep - the number of vectors to be kept at restart | ||
| 322 | |||
| 323 | Options Database Key: | ||
| 324 | . -pep_qarnoldi_restart \<keep\> - sets the restart parameter | ||
| 325 | |||
| 326 | Note: | ||
| 327 | Allowed values are in the range [0.1,0.9]. The default is 0.5. | ||
| 328 | |||
| 329 | Level: advanced | ||
| 330 | |||
| 331 | .seealso: [](ch:pep), `PEPQARNOLDI`, `PEPQArnoldiGetRestart()` | ||
| 332 | @*/ | ||
| 333 | 10 | PetscErrorCode PEPQArnoldiSetRestart(PEP pep,PetscReal keep) | |
| 334 | { | ||
| 335 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBegin; |
| 336 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
10 | PetscValidHeaderSpecific(pep,PEP_CLASSID,1); |
| 337 |
29/66✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✓ Branch 22 taken 2 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✓ Branch 26 taken 2 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✓ Branch 30 taken 2 times.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✓ Branch 33 taken 2 times.
✓ Branch 34 taken 2 times.
✗ Branch 35 not taken.
✓ Branch 36 taken 2 times.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✓ Branch 39 taken 2 times.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✓ Branch 43 taken 2 times.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✓ Branch 47 taken 2 times.
✓ Branch 48 taken 2 times.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✓ Branch 51 taken 2 times.
✓ Branch 52 taken 2 times.
✗ Branch 53 not taken.
✓ Branch 54 taken 2 times.
✗ Branch 55 not taken.
✗ Branch 56 not taken.
✓ Branch 57 taken 2 times.
✗ Branch 58 not taken.
✗ Branch 59 not taken.
✓ Branch 60 taken 2 times.
✗ Branch 61 not taken.
✗ Branch 62 not taken.
✓ Branch 63 taken 2 times.
✗ Branch 64 not taken.
✗ Branch 65 not taken.
|
10 | PetscValidLogicalCollectiveReal(pep,keep,2); |
| 338 |
8/14✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 10 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
|
10 | PetscTryMethod(pep,"PEPQArnoldiSetRestart_C",(PEP,PetscReal),(pep,keep)); |
| 339 |
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.
|
10 | PetscFunctionReturn(PETSC_SUCCESS); |
| 340 | } | ||
| 341 | |||
| 342 | 10 | static PetscErrorCode PEPQArnoldiGetRestart_QArnoldi(PEP pep,PetscReal *keep) | |
| 343 | { | ||
| 344 | 10 | PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data; | |
| 345 | |||
| 346 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBegin; |
| 347 | 10 | *keep = ctx->keep; | |
| 348 |
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.
|
10 | PetscFunctionReturn(PETSC_SUCCESS); |
| 349 | } | ||
| 350 | |||
| 351 | /*@ | ||
| 352 | PEPQArnoldiGetRestart - Gets the restart parameter used in the Q-Arnoldi method. | ||
| 353 | |||
| 354 | Not Collective | ||
| 355 | |||
| 356 | Input Parameter: | ||
| 357 | . pep - the polynomial eigensolver context | ||
| 358 | |||
| 359 | Output Parameter: | ||
| 360 | . keep - the restart parameter | ||
| 361 | |||
| 362 | Level: advanced | ||
| 363 | |||
| 364 | .seealso: [](ch:pep), `PEPQARNOLDI`, `PEPQArnoldiSetRestart()` | ||
| 365 | @*/ | ||
| 366 | 10 | PetscErrorCode PEPQArnoldiGetRestart(PEP pep,PetscReal *keep) | |
| 367 | { | ||
| 368 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBegin; |
| 369 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
10 | PetscValidHeaderSpecific(pep,PEP_CLASSID,1); |
| 370 |
2/8✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
10 | PetscAssertPointer(keep,2); |
| 371 |
9/16✓ 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 10 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
|
10 | PetscUseMethod(pep,"PEPQArnoldiGetRestart_C",(PEP,PetscReal*),(pep,keep)); |
| 372 |
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.
|
10 | PetscFunctionReturn(PETSC_SUCCESS); |
| 373 | } | ||
| 374 | |||
| 375 | 15 | static PetscErrorCode PEPQArnoldiSetLocking_QArnoldi(PEP pep,PetscBool lock) | |
| 376 | { | ||
| 377 | 15 | PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data; | |
| 378 | |||
| 379 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
15 | PetscFunctionBegin; |
| 380 | 15 | ctx->lock = lock; | |
| 381 |
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.
|
15 | PetscFunctionReturn(PETSC_SUCCESS); |
| 382 | } | ||
| 383 | |||
| 384 | /*@ | ||
| 385 | PEPQArnoldiSetLocking - Choose between locking and non-locking variants of | ||
| 386 | the Q-Arnoldi method. | ||
| 387 | |||
| 388 | Logically Collective | ||
| 389 | |||
| 390 | Input Parameters: | ||
| 391 | + pep - the polynomial eigensolver context | ||
| 392 | - lock - `PETSC_TRUE` if the locking variant must be selected | ||
| 393 | |||
| 394 | Options Database Key: | ||
| 395 | . -pep_qarnoldi_locking - sets the locking flag | ||
| 396 | |||
| 397 | Note: | ||
| 398 | The default is to lock converged eigenpairs when the method restarts. | ||
| 399 | This behavior can be changed so that all directions are kept in the | ||
| 400 | working subspace even if already converged to working accuracy (the | ||
| 401 | non-locking variant). | ||
| 402 | |||
| 403 | Level: advanced | ||
| 404 | |||
| 405 | .seealso: [](ch:pep), `PEPQARNOLDI`, `PEPQArnoldiGetLocking()` | ||
| 406 | @*/ | ||
| 407 | 15 | PetscErrorCode PEPQArnoldiSetLocking(PEP pep,PetscBool lock) | |
| 408 | { | ||
| 409 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
15 | PetscFunctionBegin; |
| 410 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
15 | PetscValidHeaderSpecific(pep,PEP_CLASSID,1); |
| 411 |
27/62✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ 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 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✓ Branch 19 taken 2 times.
✓ Branch 20 taken 2 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✓ Branch 23 taken 2 times.
✓ Branch 24 taken 2 times.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✓ Branch 27 taken 2 times.
✓ Branch 28 taken 2 times.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
✓ Branch 31 taken 2 times.
✓ Branch 32 taken 2 times.
✗ Branch 33 not taken.
✓ Branch 34 taken 2 times.
✗ Branch 35 not taken.
✗ Branch 36 not taken.
✓ Branch 37 taken 2 times.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✓ Branch 41 taken 2 times.
✗ Branch 42 not taken.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✓ Branch 45 taken 2 times.
✓ Branch 46 taken 2 times.
✗ Branch 47 not taken.
✗ Branch 48 not taken.
✓ Branch 49 taken 2 times.
✓ Branch 50 taken 2 times.
✗ Branch 51 not taken.
✓ Branch 52 taken 2 times.
✗ Branch 53 not taken.
✗ Branch 54 not taken.
✓ Branch 55 taken 2 times.
✗ Branch 56 not taken.
✗ Branch 57 not taken.
✗ Branch 58 not taken.
✓ Branch 59 taken 2 times.
✗ Branch 60 not taken.
✗ Branch 61 not taken.
|
15 | PetscValidLogicalCollectiveBool(pep,lock,2); |
| 412 |
8/14✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 10 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
|
15 | PetscTryMethod(pep,"PEPQArnoldiSetLocking_C",(PEP,PetscBool),(pep,lock)); |
| 413 |
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.
|
15 | PetscFunctionReturn(PETSC_SUCCESS); |
| 414 | } | ||
| 415 | |||
| 416 | 10 | static PetscErrorCode PEPQArnoldiGetLocking_QArnoldi(PEP pep,PetscBool *lock) | |
| 417 | { | ||
| 418 | 10 | PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data; | |
| 419 | |||
| 420 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBegin; |
| 421 | 10 | *lock = ctx->lock; | |
| 422 |
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.
|
10 | PetscFunctionReturn(PETSC_SUCCESS); |
| 423 | } | ||
| 424 | |||
| 425 | /*@ | ||
| 426 | PEPQArnoldiGetLocking - Gets the locking flag used in the Q-Arnoldi method. | ||
| 427 | |||
| 428 | Not Collective | ||
| 429 | |||
| 430 | Input Parameter: | ||
| 431 | . pep - the polynomial eigensolver context | ||
| 432 | |||
| 433 | Output Parameter: | ||
| 434 | . lock - the locking flag | ||
| 435 | |||
| 436 | Level: advanced | ||
| 437 | |||
| 438 | .seealso: [](ch:pep), `PEPQARNOLDI`, `PEPQArnoldiSetLocking()` | ||
| 439 | @*/ | ||
| 440 | 10 | PetscErrorCode PEPQArnoldiGetLocking(PEP pep,PetscBool *lock) | |
| 441 | { | ||
| 442 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBegin; |
| 443 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
10 | PetscValidHeaderSpecific(pep,PEP_CLASSID,1); |
| 444 |
2/8✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
10 | PetscAssertPointer(lock,2); |
| 445 |
9/16✓ 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 10 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
|
10 | PetscUseMethod(pep,"PEPQArnoldiGetLocking_C",(PEP,PetscBool*),(pep,lock)); |
| 446 |
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.
|
10 | PetscFunctionReturn(PETSC_SUCCESS); |
| 447 | } | ||
| 448 | |||
| 449 | 159 | static PetscErrorCode PEPSetFromOptions_QArnoldi(PEP pep,PetscOptionItems PetscOptionsObject) | |
| 450 | { | ||
| 451 | 159 | PetscBool flg,lock; | |
| 452 | 159 | PetscReal keep; | |
| 453 | |||
| 454 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
159 | PetscFunctionBegin; |
| 455 |
1/12✗ 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.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
159 | PetscOptionsHeadBegin(PetscOptionsObject,"PEP Q-Arnoldi Options"); |
| 456 | |||
| 457 |
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.
|
159 | PetscCall(PetscOptionsReal("-pep_qarnoldi_restart","Proportion of vectors kept after restart","PEPQArnoldiSetRestart",0.5,&keep,&flg)); |
| 458 |
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.
|
159 | if (flg) PetscCall(PEPQArnoldiSetRestart(pep,keep)); |
| 459 | |||
| 460 |
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.
|
159 | PetscCall(PetscOptionsBool("-pep_qarnoldi_locking","Choose between locking and non-locking variants","PEPQArnoldiSetLocking",PETSC_FALSE,&lock,&flg)); |
| 461 |
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.
|
159 | if (flg) PetscCall(PEPQArnoldiSetLocking(pep,lock)); |
| 462 | |||
| 463 |
1/14✗ Branch 0 not taken.
✓ Branch 1 taken 2 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.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
159 | PetscOptionsHeadEnd(); |
| 464 |
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.
|
31 | PetscFunctionReturn(PETSC_SUCCESS); |
| 465 | } | ||
| 466 | |||
| 467 | ✗ | static PetscErrorCode PEPView_QArnoldi(PEP pep,PetscViewer viewer) | |
| 468 | { | ||
| 469 | ✗ | PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data; | |
| 470 | ✗ | PetscBool isascii; | |
| 471 | |||
| 472 | ✗ | PetscFunctionBegin; | |
| 473 | ✗ | PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii)); | |
| 474 | ✗ | if (isascii) { | |
| 475 | ✗ | PetscCall(PetscViewerASCIIPrintf(viewer," %d%% of basis vectors kept after restart\n",(int)(100*ctx->keep))); | |
| 476 | ✗ | PetscCall(PetscViewerASCIIPrintf(viewer," using the %slocking variant\n",ctx->lock?"":"non-")); | |
| 477 | } | ||
| 478 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
| 479 | } | ||
| 480 | |||
| 481 | 169 | static PetscErrorCode PEPDestroy_QArnoldi(PEP pep) | |
| 482 | { | ||
| 483 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
169 | PetscFunctionBegin; |
| 484 |
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.
|
169 | PetscCall(PetscFree(pep->data)); |
| 485 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetRestart_C",NULL)); |
| 486 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetRestart_C",NULL)); |
| 487 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetLocking_C",NULL)); |
| 488 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetLocking_C",NULL)); |
| 489 |
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.
|
33 | PetscFunctionReturn(PETSC_SUCCESS); |
| 490 | } | ||
| 491 | |||
| 492 | /*MC | ||
| 493 | PEPQARNOLDI - PEPQARNOLDI = "qarnoldi" - The Quadratic Arnoldi method (Q-Arnoldi) | ||
| 494 | for quadratic eigenvalue problems. | ||
| 495 | |||
| 496 | Notes: | ||
| 497 | This solver is available for quadratic eigenproblems only. | ||
| 498 | |||
| 499 | It implements the Q-Arnoldi method {cite:p}`Mee09`, which is not guaranteed | ||
| 500 | to be numerically stable. Users should prefer `PEPTOAR`, which | ||
| 501 | provides a similar algorithm with guaranteed numerical stability. | ||
| 502 | |||
| 503 | Level: beginner | ||
| 504 | |||
| 505 | .seealso: [](ch:pep), `PEP`, `PEPType`, `PEPSetType()` | ||
| 506 | M*/ | ||
| 507 | 169 | SLEPC_EXTERN PetscErrorCode PEPCreate_QArnoldi(PEP pep) | |
| 508 | { | ||
| 509 | 169 | PEP_QARNOLDI *ctx; | |
| 510 | |||
| 511 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
169 | PetscFunctionBegin; |
| 512 |
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.
|
169 | PetscCall(PetscNew(&ctx)); |
| 513 | 169 | pep->data = (void*)ctx; | |
| 514 | |||
| 515 | 169 | pep->lineariz = PETSC_TRUE; | |
| 516 | 169 | ctx->lock = PETSC_TRUE; | |
| 517 | |||
| 518 | 169 | pep->ops->solve = PEPSolve_QArnoldi; | |
| 519 | 169 | pep->ops->setup = PEPSetUp_QArnoldi; | |
| 520 | 169 | pep->ops->setfromoptions = PEPSetFromOptions_QArnoldi; | |
| 521 | 169 | pep->ops->destroy = PEPDestroy_QArnoldi; | |
| 522 | 169 | pep->ops->view = PEPView_QArnoldi; | |
| 523 | 169 | pep->ops->backtransform = PEPBackTransform_Default; | |
| 524 | 169 | pep->ops->computevectors = PEPComputeVectors_Default; | |
| 525 | 169 | pep->ops->extractvectors = PEPExtractVectors_QArnoldi; | |
| 526 | 169 | pep->ops->setdefaultst = PEPSetDefaultST_Transform; | |
| 527 | |||
| 528 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetRestart_C",PEPQArnoldiSetRestart_QArnoldi)); |
| 529 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetRestart_C",PEPQArnoldiGetRestart_QArnoldi)); |
| 530 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetLocking_C",PEPQArnoldiSetLocking_QArnoldi)); |
| 531 |
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.
|
169 | PetscCall(PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetLocking_C",PEPQArnoldiGetLocking_QArnoldi)); |
| 532 |
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.
|
33 | PetscFunctionReturn(PETSC_SUCCESS); |
| 533 | } | ||
| 534 |