| 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 | Two-sided variant of the NEPSLP solver. | ||
| 12 | */ | ||
| 13 | |||
| 14 | #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/ | ||
| 15 | #include "slp.h" | ||
| 16 | |||
| 17 | typedef struct _n_nep_def_ctx *NEP_NEDEF_CTX; | ||
| 18 | |||
| 19 | struct _n_nep_def_ctx { | ||
| 20 | PetscInt n; | ||
| 21 | PetscBool ref; | ||
| 22 | PetscScalar *eig; | ||
| 23 | BV V,W; | ||
| 24 | }; | ||
| 25 | |||
| 26 | typedef struct { /* context for two-sided solver */ | ||
| 27 | Mat Ft; | ||
| 28 | Mat Jt; | ||
| 29 | Vec w; | ||
| 30 | } NEP_SLPTS_MATSHELL; | ||
| 31 | |||
| 32 | typedef struct { /* context for non-equivalence deflation */ | ||
| 33 | NEP_NEDEF_CTX defctx; | ||
| 34 | Mat F; | ||
| 35 | Mat P; | ||
| 36 | Mat J; | ||
| 37 | KSP ksp; | ||
| 38 | PetscBool isJ; | ||
| 39 | PetscScalar lambda; | ||
| 40 | Vec w[2]; | ||
| 41 | } NEP_NEDEF_MATSHELL; | ||
| 42 | |||
| 43 | 4760 | static PetscErrorCode MatMult_SLPTS_Right(Mat M,Vec x,Vec y) | |
| 44 | { | ||
| 45 | 4760 | NEP_SLPTS_MATSHELL *ctx; | |
| 46 | |||
| 47 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
4760 | PetscFunctionBegin; |
| 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.
|
4760 | PetscCall(MatShellGetContext(M,&ctx)); |
| 49 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
4760 | PetscCall(MatMult(ctx->Jt,x,ctx->w)); |
| 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.
|
4760 | PetscCall(MatSolve(ctx->Ft,ctx->w,y)); |
| 51 |
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.
|
952 | PetscFunctionReturn(PETSC_SUCCESS); |
| 52 | } | ||
| 53 | |||
| 54 | 4960 | static PetscErrorCode MatMult_SLPTS_Left(Mat M,Vec x,Vec y) | |
| 55 | { | ||
| 56 | 4960 | NEP_SLPTS_MATSHELL *ctx; | |
| 57 | |||
| 58 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
4960 | PetscFunctionBegin; |
| 59 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
4960 | PetscCall(MatShellGetContext(M,&ctx)); |
| 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.
|
4960 | PetscCall(MatMultTranspose(ctx->Jt,x,ctx->w)); |
| 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.
|
4960 | PetscCall(MatSolveTranspose(ctx->Ft,ctx->w,y)); |
| 62 |
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.
|
992 | PetscFunctionReturn(PETSC_SUCCESS); |
| 63 | } | ||
| 64 | |||
| 65 | 100 | static PetscErrorCode MatDestroy_SLPTS(Mat M) | |
| 66 | { | ||
| 67 | 100 | NEP_SLPTS_MATSHELL *ctx; | |
| 68 | |||
| 69 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
100 | PetscFunctionBegin; |
| 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.
|
100 | PetscCall(MatShellGetContext(M,&ctx)); |
| 71 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(VecDestroy(&ctx->w)); |
| 72 |
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.
|
100 | PetscCall(PetscFree(ctx)); |
| 73 |
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.
|
20 | PetscFunctionReturn(PETSC_SUCCESS); |
| 74 | } | ||
| 75 | |||
| 76 | #if defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP) | ||
| 77 | 84 | static PetscErrorCode MatCreateVecs_SLPTS(Mat M,Vec *left,Vec *right) | |
| 78 | { | ||
| 79 | 84 | NEP_SLPTS_MATSHELL *ctx; | |
| 80 | |||
| 81 | 84 | PetscFunctionBegin; | |
| 82 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 4 times.
|
84 | PetscCall(MatShellGetContext(M,&ctx)); |
| 83 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
84 | if (right) PetscCall(VecDuplicate(ctx->w,right)); |
| 84 |
2/4✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 4 times.
|
84 | if (left) PetscCall(VecDuplicate(ctx->w,left)); |
| 85 | PetscFunctionReturn(PETSC_SUCCESS); | ||
| 86 | } | ||
| 87 | #endif | ||
| 88 | |||
| 89 | 100 | static PetscErrorCode NEPSLPSetUpEPSMat(NEP nep,Mat F,Mat J,PetscBool left,Mat *M) | |
| 90 | { | ||
| 91 | 100 | Mat Mshell; | |
| 92 | 100 | PetscInt nloc,mloc; | |
| 93 | 100 | NEP_SLPTS_MATSHELL *shellctx; | |
| 94 | |||
| 95 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
100 | PetscFunctionBegin; |
| 96 | /* Create mat shell */ | ||
| 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.
|
100 | PetscCall(PetscNew(&shellctx)); |
| 98 | 100 | shellctx->Ft = F; | |
| 99 | 100 | shellctx->Jt = J; | |
| 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.
|
100 | PetscCall(MatGetLocalSize(nep->function,&mloc,&nloc)); |
| 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.
|
100 | PetscCall(MatCreateShell(PetscObjectComm((PetscObject)nep),nloc,mloc,PETSC_DETERMINE,PETSC_DETERMINE,shellctx,&Mshell)); |
| 102 |
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.
|
100 | if (left) PetscCall(MatShellSetOperation(Mshell,MATOP_MULT,(PetscErrorCodeFn*)MatMult_SLPTS_Left)); |
| 103 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | else PetscCall(MatShellSetOperation(Mshell,MATOP_MULT,(PetscErrorCodeFn*)MatMult_SLPTS_Right)); |
| 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.
|
100 | PetscCall(MatShellSetOperation(Mshell,MATOP_DESTROY,(PetscErrorCodeFn*)MatDestroy_SLPTS)); |
| 105 | #if defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP) | ||
| 106 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 4 times.
|
40 | PetscCall(MatShellSetOperation(Mshell,MATOP_CREATE_VECS,(PetscErrorCodeFn*)MatCreateVecs_SLPTS)); |
| 107 | #endif | ||
| 108 | 100 | *M = Mshell; | |
| 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.
|
100 | PetscCall(MatCreateVecs(nep->function,&shellctx->w,NULL)); |
| 110 |
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.
|
20 | PetscFunctionReturn(PETSC_SUCCESS); |
| 111 | } | ||
| 112 | |||
| 113 | /* Functions for deflation */ | ||
| 114 | 50 | static PetscErrorCode NEPDeflationNEDestroy(NEP_NEDEF_CTX defctx) | |
| 115 | { | ||
| 116 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
50 | PetscFunctionBegin; |
| 117 |
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.
|
50 | if (!defctx) PetscFunctionReturn(PETSC_SUCCESS); |
| 118 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
50 | PetscCall(PetscFree(defctx->eig)); |
| 119 |
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.
|
50 | PetscCall(PetscFree(defctx)); |
| 120 |
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); |
| 121 | } | ||
| 122 | |||
| 123 | 50 | static PetscErrorCode NEPDeflationNECreate(NEP nep,BV V,BV W,PetscInt sz,NEP_NEDEF_CTX *defctx) | |
| 124 | { | ||
| 125 | 50 | NEP_NEDEF_CTX op; | |
| 126 | |||
| 127 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
50 | PetscFunctionBegin; |
| 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.
|
50 | PetscCall(PetscNew(&op)); |
| 129 | 50 | *defctx = op; | |
| 130 | 50 | op->n = 0; | |
| 131 | 50 | op->ref = PETSC_FALSE; | |
| 132 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(PetscCalloc1(sz,&op->eig)); |
| 133 |
2/4✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
|
50 | PetscCall(PetscObjectStateIncrease((PetscObject)V)); |
| 134 |
2/4✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
|
50 | PetscCall(PetscObjectStateIncrease((PetscObject)W)); |
| 135 | 50 | op->V = V; | |
| 136 | 50 | op->W = W; | |
| 137 |
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.
|
50 | PetscFunctionReturn(PETSC_SUCCESS); |
| 138 | } | ||
| 139 | |||
| 140 | 570 | static PetscErrorCode NEPDeflationNEComputeFunction(NEP nep,Mat M,PetscScalar lambda) | |
| 141 | { | ||
| 142 | 570 | NEP_NEDEF_MATSHELL *matctx; | |
| 143 | |||
| 144 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
570 | PetscFunctionBegin; |
| 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.
|
570 | PetscCall(MatShellGetContext(M,&matctx)); |
| 146 |
8/14✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 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.
|
570 | if (lambda==matctx->lambda) PetscFunctionReturn(PETSC_SUCCESS); |
| 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.
|
410 | PetscCall(NEPComputeFunction(nep,lambda,matctx->F,matctx->P)); |
| 148 |
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.
|
410 | if (matctx->isJ) PetscCall(NEPComputeJacobian(nep,lambda,matctx->J)); |
| 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.
|
410 | if (matctx->ksp) PetscCall(NEP_KSPSetOperators(matctx->ksp,matctx->F,matctx->P)); |
| 150 | 410 | matctx->lambda = lambda; | |
| 151 |
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.
|
410 | PetscFunctionReturn(PETSC_SUCCESS); |
| 152 | } | ||
| 153 | |||
| 154 | 5050 | static PetscErrorCode MatMult_NEPDeflationNE(Mat M,Vec x,Vec r) | |
| 155 | { | ||
| 156 | 5050 | Vec t,tt; | |
| 157 | 5050 | PetscScalar *h,*alpha,lambda,*eig; | |
| 158 | 5050 | PetscInt i,k; | |
| 159 | 5050 | NEP_NEDEF_MATSHELL *matctx; | |
| 160 | |||
| 161 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
5050 | PetscFunctionBegin; |
| 162 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5050 | PetscCall(MatShellGetContext(M,&matctx)); |
| 163 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
5050 | if (matctx->defctx->n && !matctx->defctx->ref) { |
| 164 | 890 | k = matctx->defctx->n; | |
| 165 | 890 | lambda = matctx->lambda; | |
| 166 | 890 | eig = matctx->defctx->eig; | |
| 167 | 890 | t = matctx->w[0]; | |
| 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.
|
890 | PetscCall(VecCopy(x,t)); |
| 169 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
890 | PetscCall(PetscMalloc2(k,&h,k,&alpha)); |
| 170 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1960 | for (i=0;i<k;i++) alpha[i] = (lambda-eig[i]-1.0)/(lambda-eig[i]); |
| 171 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
890 | PetscCall(BVDotVec(matctx->defctx->V,t,h)); |
| 172 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1960 | for (i=0;i<k;i++) h[i] *= alpha[i]; |
| 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.
|
890 | PetscCall(BVMultVec(matctx->defctx->W,-1.0,1.0,t,h)); |
| 174 |
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.
|
890 | PetscCall(MatMult(matctx->isJ?matctx->J:matctx->F,t,r)); |
| 175 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
890 | if (matctx->isJ) { |
| 176 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1760 | for (i=0;i<k;i++) h[i] *= (1.0/((lambda-eig[i])*(lambda-eig[i])))/alpha[i]; |
| 177 | 800 | tt = matctx->w[1]; | |
| 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.
|
800 | PetscCall(BVMultVec(matctx->defctx->W,-1.0,0.0,tt,h)); |
| 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.
|
800 | PetscCall(MatMult(matctx->F,tt,t)); |
| 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.
|
800 | PetscCall(VecAXPY(r,1.0,t)); |
| 181 | } | ||
| 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.
|
890 | PetscCall(PetscFree2(h,alpha)); |
| 183 |
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.
|
4160 | } else PetscCall(MatMult(matctx->isJ?matctx->J:matctx->F,x,r)); |
| 184 |
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.
|
1010 | PetscFunctionReturn(PETSC_SUCCESS); |
| 185 | } | ||
| 186 | |||
| 187 | 5250 | static PetscErrorCode MatMultTranspose_NEPDeflationNE(Mat M,Vec x,Vec r) | |
| 188 | { | ||
| 189 | 5250 | Vec t,tt; | |
| 190 | 5250 | PetscScalar *h,*alphaC,lambda,*eig; | |
| 191 | 5250 | PetscInt i,k; | |
| 192 | 5250 | NEP_NEDEF_MATSHELL *matctx; | |
| 193 | |||
| 194 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
5250 | PetscFunctionBegin; |
| 195 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5250 | PetscCall(MatShellGetContext(M,&matctx)); |
| 196 | 5250 | t = matctx->w[0]; | |
| 197 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5250 | PetscCall(VecCopy(x,t)); |
| 198 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
5250 | if (matctx->defctx->n && !matctx->defctx->ref) { |
| 199 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
890 | PetscCall(VecConjugate(t)); |
| 200 | 890 | k = matctx->defctx->n; | |
| 201 | 890 | lambda = matctx->lambda; | |
| 202 | 890 | eig = matctx->defctx->eig; | |
| 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.
|
890 | PetscCall(PetscMalloc2(k,&h,k,&alphaC)); |
| 204 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1960 | for (i=0;i<k;i++) alphaC[i] = PetscConj((lambda-eig[i]-1.0)/(lambda-eig[i])); |
| 205 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
890 | PetscCall(BVDotVec(matctx->defctx->W,t,h)); |
| 206 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1960 | for (i=0;i<k;i++) h[i] *= alphaC[i]; |
| 207 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
890 | PetscCall(BVMultVec(matctx->defctx->V,-1.0,1.0,t,h)); |
| 208 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
890 | PetscCall(VecConjugate(t)); |
| 209 |
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.
|
890 | PetscCall(MatMultTranspose(matctx->isJ?matctx->J:matctx->F,t,r)); |
| 210 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
890 | if (matctx->isJ) { |
| 211 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1760 | for (i=0;i<k;i++) h[i] *= PetscConj(1.0/((lambda-eig[i])*(lambda-eig[i])))/alphaC[i]; |
| 212 | 800 | tt = matctx->w[1]; | |
| 213 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(BVMultVec(matctx->defctx->V,-1.0,0.0,tt,h)); |
| 214 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(VecConjugate(tt)); |
| 215 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(MatMultTranspose(matctx->F,tt,t)); |
| 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.
|
800 | PetscCall(VecAXPY(r,1.0,t)); |
| 217 | } | ||
| 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.
|
890 | PetscCall(PetscFree2(h,alphaC)); |
| 219 |
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.
|
4360 | } else PetscCall(MatMultTranspose(matctx->isJ?matctx->J:matctx->F,t,r)); |
| 220 |
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.
|
1050 | PetscFunctionReturn(PETSC_SUCCESS); |
| 221 | } | ||
| 222 | |||
| 223 | 4760 | static PetscErrorCode MatSolve_NEPDeflationNE(Mat M,Vec b,Vec x) | |
| 224 | { | ||
| 225 | 4760 | PetscScalar *h,*alpha,lambda,*eig; | |
| 226 | 4760 | PetscInt i,k; | |
| 227 | 4760 | NEP_NEDEF_MATSHELL *matctx; | |
| 228 | |||
| 229 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
4760 | PetscFunctionBegin; |
| 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.
|
4760 | PetscCall(MatShellGetContext(M,&matctx)); |
| 231 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
4760 | if (!matctx->ksp) { |
| 232 | ✗ | PetscCall(VecCopy(b,x)); | |
| 233 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
| 234 | } | ||
| 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.
|
4760 | PetscCall(KSPSolve(matctx->ksp,b,x)); |
| 236 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
4760 | if (matctx->defctx->n && !matctx->defctx->ref) { |
| 237 | 800 | k = matctx->defctx->n; | |
| 238 | 800 | lambda = matctx->lambda; | |
| 239 | 800 | eig = matctx->defctx->eig; | |
| 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.
|
800 | PetscCall(PetscMalloc2(k,&h,k,&alpha)); |
| 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.
|
800 | PetscCall(BVDotVec(matctx->defctx->V,x,h)); |
| 242 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1760 | for (i=0;i<k;i++) alpha[i] = (lambda-eig[i]-1.0)/(lambda-eig[i]); |
| 243 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1760 | for (i=0;i<k;i++) h[i] *= alpha[i]/(1.0-alpha[i]); |
| 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.
|
800 | PetscCall(BVMultVec(matctx->defctx->W,1.0,1.0,x,h)); |
| 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.
|
800 | PetscCall(PetscFree2(h,alpha)); |
| 246 | } | ||
| 247 |
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.
|
952 | PetscFunctionReturn(PETSC_SUCCESS); |
| 248 | } | ||
| 249 | |||
| 250 | 4960 | static PetscErrorCode MatSolveTranspose_NEPDeflationNE(Mat M,Vec b,Vec x) | |
| 251 | { | ||
| 252 | 4960 | PetscScalar *h,*alphaC,lambda,*eig; | |
| 253 | 4960 | PetscInt i,k; | |
| 254 | 4960 | NEP_NEDEF_MATSHELL *matctx; | |
| 255 | |||
| 256 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
4960 | PetscFunctionBegin; |
| 257 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
4960 | PetscCall(MatShellGetContext(M,&matctx)); |
| 258 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
4960 | if (!matctx->ksp) { |
| 259 | ✗ | PetscCall(VecCopy(b,x)); | |
| 260 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
| 261 | } | ||
| 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.
|
4960 | PetscCall(KSPSolveTranspose(matctx->ksp,b,x)); |
| 263 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
4960 | if (matctx->defctx->n && !matctx->defctx->ref) { |
| 264 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(VecConjugate(x)); |
| 265 | 800 | k = matctx->defctx->n; | |
| 266 | 800 | lambda = matctx->lambda; | |
| 267 | 800 | eig = matctx->defctx->eig; | |
| 268 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(PetscMalloc2(k,&h,k,&alphaC)); |
| 269 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(BVDotVec(matctx->defctx->W,x,h)); |
| 270 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1760 | for (i=0;i<k;i++) alphaC[i] = PetscConj((lambda-eig[i]-1.0)/(lambda-eig[i])); |
| 271 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1760 | for (i=0;i<k;i++) h[i] *= alphaC[i]/(1.0-alphaC[i]); |
| 272 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(BVMultVec(matctx->defctx->V,1.0,1.0,x,h)); |
| 273 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
800 | PetscCall(PetscFree2(h,alphaC)); |
| 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.
|
800 | PetscCall(VecConjugate(x)); |
| 275 | } | ||
| 276 |
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.
|
992 | PetscFunctionReturn(PETSC_SUCCESS); |
| 277 | } | ||
| 278 | |||
| 279 | 100 | static PetscErrorCode MatDestroy_NEPDeflationNE(Mat M) | |
| 280 | { | ||
| 281 | 100 | NEP_NEDEF_MATSHELL *matctx; | |
| 282 | |||
| 283 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
100 | PetscFunctionBegin; |
| 284 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatShellGetContext(M,&matctx)); |
| 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.
|
100 | PetscCall(VecDestroy(&matctx->w[0])); |
| 286 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(VecDestroy(&matctx->w[1])); |
| 287 |
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.
|
100 | PetscCall(PetscFree(matctx)); |
| 288 |
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.
|
20 | PetscFunctionReturn(PETSC_SUCCESS); |
| 289 | } | ||
| 290 | |||
| 291 | ✗ | static PetscErrorCode MatCreateVecs_NEPDeflationNE(Mat M,Vec *right,Vec *left) | |
| 292 | { | ||
| 293 | ✗ | NEP_NEDEF_MATSHELL *matctx; | |
| 294 | |||
| 295 | ✗ | PetscFunctionBegin; | |
| 296 | ✗ | PetscCall(MatShellGetContext(M,&matctx)); | |
| 297 | ✗ | PetscCall(MatCreateVecs(matctx->F,right,left)); | |
| 298 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
| 299 | } | ||
| 300 | |||
| 301 | 100 | static PetscErrorCode NEPDeflationNEFunctionCreate(NEP_NEDEF_CTX defctx,NEP nep,Mat F,Mat P,Mat J,KSP ksp,PetscBool isJ,Mat *Mshell) | |
| 302 | { | ||
| 303 | 100 | NEP_NEDEF_MATSHELL *matctx; | |
| 304 | 100 | PetscInt nloc,mloc; | |
| 305 | |||
| 306 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
100 | PetscFunctionBegin; |
| 307 | /* Create mat shell */ | ||
| 308 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(PetscNew(&matctx)); |
| 309 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatGetLocalSize(nep->function,&mloc,&nloc)); |
| 310 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatCreateShell(PetscObjectComm((PetscObject)nep),nloc,mloc,PETSC_DETERMINE,PETSC_DETERMINE,matctx,Mshell)); |
| 311 | 100 | matctx->F = F; | |
| 312 | 100 | matctx->P = P; | |
| 313 | 100 | matctx->J = J; | |
| 314 | 100 | matctx->isJ = isJ; | |
| 315 | 100 | matctx->ksp = ksp; | |
| 316 | 100 | matctx->defctx = defctx; | |
| 317 | 100 | matctx->lambda = PETSC_MAX_REAL; | |
| 318 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatCreateVecs(F,&matctx->w[0],NULL)); |
| 319 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(VecDuplicate(matctx->w[0],&matctx->w[1])); |
| 320 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatShellSetOperation(*Mshell,MATOP_MULT,(PetscErrorCodeFn*)MatMult_NEPDeflationNE)); |
| 321 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatShellSetOperation(*Mshell,MATOP_MULT_TRANSPOSE,(PetscErrorCodeFn*)MatMultTranspose_NEPDeflationNE)); |
| 322 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatShellSetOperation(*Mshell,MATOP_SOLVE,(PetscErrorCodeFn*)MatSolve_NEPDeflationNE)); |
| 323 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatShellSetOperation(*Mshell,MATOP_SOLVE_TRANSPOSE,(PetscErrorCodeFn*)MatSolveTranspose_NEPDeflationNE)); |
| 324 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatShellSetOperation(*Mshell,MATOP_DESTROY,(PetscErrorCodeFn*)MatDestroy_NEPDeflationNE)); |
| 325 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
100 | PetscCall(MatShellSetOperation(*Mshell,MATOP_CREATE_VECS,(PetscErrorCodeFn*)MatCreateVecs_NEPDeflationNE)); |
| 326 |
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.
|
20 | PetscFunctionReturn(PETSC_SUCCESS); |
| 327 | } | ||
| 328 | |||
| 329 | 130 | static PetscErrorCode NEPDeflationNERecoverEigenvectors(NEP_NEDEF_CTX defctx,Vec u,Vec w,PetscScalar lambda) | |
| 330 | { | ||
| 331 | 130 | PetscScalar *h,*alpha,*eig; | |
| 332 | 130 | PetscInt i,k; | |
| 333 | |||
| 334 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
130 | PetscFunctionBegin; |
| 335 |
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.
|
130 | if (w) PetscCall(VecConjugate(w)); |
| 336 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
130 | if (defctx->n && !defctx->ref) { |
| 337 | 40 | eig = defctx->eig; | |
| 338 | 40 | k = defctx->n; | |
| 339 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(PetscMalloc2(k,&h,k,&alpha)); |
| 340 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
90 | for (i=0;i<k;i++) alpha[i] = (lambda-eig[i]-1.0)/(lambda-eig[i]); |
| 341 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(BVDotVec(defctx->V,u,h)); |
| 342 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
90 | for (i=0;i<k;i++) h[i] *= alpha[i]; |
| 343 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(BVMultVec(defctx->W,-1.0,1.0,u,h)); |
| 344 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(VecNormalize(u,NULL)); |
| 345 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
40 | if (w) { |
| 346 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(BVDotVec(defctx->W,w,h)); |
| 347 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
90 | for (i=0;i<k;i++) alpha[i] = PetscConj((lambda-eig[i]-1.0)/(lambda-eig[i])); |
| 348 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
90 | for (i=0;i<k;i++) h[i] *= alpha[i]; |
| 349 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(BVMultVec(defctx->V,-1.0,1.0,w,h)); |
| 350 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(VecNormalize(w,NULL)); |
| 351 | } | ||
| 352 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(PetscFree2(h,alpha)); |
| 353 | } | ||
| 354 |
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.
|
26 | PetscFunctionReturn(PETSC_SUCCESS); |
| 355 | } | ||
| 356 | |||
| 357 | 90 | static PetscErrorCode NEPDeflationNELocking(NEP_NEDEF_CTX defctx,Vec u,Vec w,PetscScalar lambda) | |
| 358 | { | ||
| 359 | 90 | PetscInt n; | |
| 360 | |||
| 361 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
90 | PetscFunctionBegin; |
| 362 | 90 | n = defctx->n++; | |
| 363 | 90 | defctx->eig[n] = lambda; | |
| 364 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(BVInsertVec(defctx->V,n,u)); |
| 365 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(BVInsertVec(defctx->W,n,w)); |
| 366 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(BVSetActiveColumns(defctx->V,0,defctx->n)); |
| 367 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(BVSetActiveColumns(defctx->W,0,defctx->n)); |
| 368 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(BVBiorthonormalizeColumn(defctx->V,defctx->W,n,NULL)); |
| 369 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
18 | PetscFunctionReturn(PETSC_SUCCESS); |
| 370 | } | ||
| 371 | |||
| 372 | 130 | static PetscErrorCode NEPDeflationNESetRefine(NEP_NEDEF_CTX defctx,PetscBool ref) | |
| 373 | { | ||
| 374 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
130 | PetscFunctionBegin; |
| 375 | 130 | defctx->ref = ref; | |
| 376 |
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.
|
130 | PetscFunctionReturn(PETSC_SUCCESS); |
| 377 | } | ||
| 378 | |||
| 379 | 50 | PetscErrorCode NEPSolve_SLP_Twosided(NEP nep) | |
| 380 | { | ||
| 381 | 50 | NEP_SLP *ctx = (NEP_SLP*)nep->data; | |
| 382 | 50 | Mat mF,mJ,M,Mt; | |
| 383 | 50 | Vec u,r,t,w; | |
| 384 | 50 | BV X,Y; | |
| 385 | 50 | PetscScalar sigma,lambda,mu,im=0.0,mu2,im2; | |
| 386 | 50 | PetscReal resnorm,resl; | |
| 387 | 50 | PetscInt nconv,nconv2,i; | |
| 388 | 50 | PetscBool skip=PETSC_FALSE,lock=PETSC_FALSE; | |
| 389 | 50 | NEP_NEDEF_CTX defctx=NULL; /* Extended operator for deflation */ | |
| 390 | |||
| 391 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
50 | PetscFunctionBegin; |
| 392 | /* get initial approximation of eigenvalue and eigenvector */ | ||
| 393 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(NEPGetDefaultShift(nep,&sigma)); |
| 394 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
50 | if (!nep->nini) PetscCall(BVSetRandomColumn(nep->V,0)); |
| 395 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVSetRandomColumn(nep->W,0)); |
| 396 | 50 | lambda = sigma; | |
| 397 |
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.
|
50 | if (!ctx->ksp) PetscCall(NEPSLPGetKSP(nep,&ctx->ksp)); |
| 398 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVDuplicate(nep->V,&X)); |
| 399 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVDuplicate(nep->W,&Y)); |
| 400 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(NEPDeflationNECreate(nep,X,Y,nep->nev,&defctx)); |
| 401 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVGetColumn(nep->V,0,&t)); |
| 402 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(VecDuplicate(t,&u)); |
| 403 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(VecDuplicate(t,&w)); |
| 404 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVRestoreColumn(nep->V,0,&t)); |
| 405 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVCopyVec(nep->V,0,u)); |
| 406 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVCopyVec(nep->W,0,w)); |
| 407 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(VecDuplicate(u,&r)); |
| 408 |
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.
|
50 | PetscCall(NEPDeflationNEFunctionCreate(defctx,nep,nep->function,nep->function_pre?nep->function_pre:nep->function,NULL,ctx->ksp,PETSC_FALSE,&mF)); |
| 409 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(NEPDeflationNEFunctionCreate(defctx,nep,nep->function,nep->function,nep->jacobian,NULL,PETSC_TRUE,&mJ)); |
| 410 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(NEPSLPSetUpEPSMat(nep,mF,mJ,PETSC_FALSE,&M)); |
| 411 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(NEPSLPSetUpEPSMat(nep,mF,mJ,PETSC_TRUE,&Mt)); |
| 412 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(EPSSetOperators(ctx->eps,M,NULL)); |
| 413 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(MatDestroy(&M)); |
| 414 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(EPSSetOperators(ctx->epsts,Mt,NULL)); |
| 415 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(MatDestroy(&Mt)); |
| 416 | |||
| 417 | /* Restart loop */ | ||
| 418 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
300 | while (nep->reason == NEP_CONVERGED_ITERATING) { |
| 419 | 250 | nep->its++; | |
| 420 | |||
| 421 | /* form residual, r = T(lambda)*u (used in convergence test only) */ | ||
| 422 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
250 | PetscCall(NEPDeflationNEComputeFunction(nep,mF,lambda)); |
| 423 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
250 | PetscCall(MatMultTranspose(mF,w,r)); |
| 424 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
250 | PetscCall(VecNorm(r,NORM_2,&resl)); |
| 425 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
250 | PetscCall(MatMult(mF,u,r)); |
| 426 | |||
| 427 | /* convergence test */ | ||
| 428 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
250 | PetscCall(VecNorm(r,NORM_2,&resnorm)); |
| 429 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
250 | resnorm = PetscMax(resnorm,resl); |
| 430 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
250 | PetscCall((*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx)); |
| 431 | 250 | nep->eigr[nep->nconv] = lambda; | |
| 432 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
250 | if (nep->errest[nep->nconv]<=nep->tol || nep->errest[nep->nconv]<=ctx->deftol) { |
| 433 |
4/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 10 times.
✓ Branch 5 taken 10 times.
|
90 | if (nep->errest[nep->nconv]<=ctx->deftol && !defctx->ref && nep->nconv) { |
| 434 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(NEPDeflationNERecoverEigenvectors(defctx,u,w,lambda)); |
| 435 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(VecConjugate(w)); |
| 436 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(NEPDeflationNESetRefine(defctx,PETSC_TRUE)); |
| 437 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(MatMultTranspose(mF,w,r)); |
| 438 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(VecNorm(r,NORM_2,&resl)); |
| 439 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(MatMult(mF,u,r)); |
| 440 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(VecNorm(r,NORM_2,&resnorm)); |
| 441 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 9 times.
|
40 | resnorm = PetscMax(resnorm,resl); |
| 442 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall((*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx)); |
| 443 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
40 | if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE; |
| 444 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
50 | } else if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE; |
| 445 | } | ||
| 446 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
160 | if (lock) { |
| 447 | 90 | lock = PETSC_FALSE; | |
| 448 | 90 | skip = PETSC_TRUE; | |
| 449 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(NEPDeflationNERecoverEigenvectors(defctx,u,w,lambda)); |
| 450 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(NEPDeflationNELocking(defctx,u,w,lambda)); |
| 451 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(NEPDeflationNESetRefine(defctx,PETSC_FALSE)); |
| 452 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(BVInsertVec(nep->V,nep->nconv,u)); |
| 453 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(BVInsertVec(nep->W,nep->nconv,w)); |
| 454 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
90 | PetscCall(VecConjugate(w)); |
| 455 | 90 | nep->nconv = nep->nconv + 1; | |
| 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.
|
250 | PetscCall((*nep->stopping)(nep,nep->its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx)); |
| 458 |
11/12✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
✓ Branch 4 taken 10 times.
✓ Branch 5 taken 8 times.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 8 times.
✓ Branch 8 taken 2 times.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
250 | if (!skip || nep->reason>0) PetscCall(NEPMonitor(nep,nep->its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1)); |
| 459 | |||
| 460 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
250 | if (nep->reason == NEP_CONVERGED_ITERATING) { |
| 461 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
200 | if (!skip) { |
| 462 | /* evaluate T(lambda) and T'(lambda) */ | ||
| 463 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(NEPDeflationNEComputeFunction(nep,mF,lambda)); |
| 464 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(NEPDeflationNEComputeFunction(nep,mJ,lambda)); |
| 465 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSSetInitialSpace(ctx->eps,1,&u)); |
| 466 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSSetInitialSpace(ctx->epsts,1,&w)); |
| 467 | |||
| 468 | /* compute new eigenvalue correction mu and eigenvector approximation u */ | ||
| 469 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSSolve(ctx->eps)); |
| 470 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSSolve(ctx->epsts)); |
| 471 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSGetConverged(ctx->eps,&nconv)); |
| 472 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSGetConverged(ctx->epsts,&nconv2)); |
| 473 |
2/4✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
160 | if (!nconv||!nconv2) { |
| 474 | ✗ | PetscCall(PetscInfo(nep,"iter=%" PetscInt_FMT ", inner iteration failed, stopping solve\n",nep->its)); | |
| 475 | ✗ | nep->reason = NEP_DIVERGED_LINEAR_SOLVE; | |
| 476 | ✗ | break; | |
| 477 | } | ||
| 478 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSGetEigenpair(ctx->eps,0,&mu,&im,u,NULL)); |
| 479 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
160 | for (i=0;i<nconv2;i++) { |
| 480 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
160 | PetscCall(EPSGetEigenpair(ctx->epsts,i,&mu2,&im2,w,NULL)); |
| 481 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
160 | if (SlepcAbsEigenvalue(mu-mu2,im-im2)/SlepcAbsEigenvalue(mu,im)<nep->tol*1000) break; |
| 482 | } | ||
| 483 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
160 | if (i==nconv2) { |
| 484 | ✗ | PetscCall(PetscInfo(nep,"iter=%" PetscInt_FMT ", inner iteration failed, stopping solve\n",nep->its)); | |
| 485 | ✗ | nep->reason = NEP_DIVERGED_LINEAR_SOLVE; | |
| 486 | ✗ | break; | |
| 487 | } | ||
| 488 | |||
| 489 | 160 | mu = 1.0/mu; | |
| 490 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
160 | PetscCheck(PetscAbsScalar(im)<PETSC_MACHINE_EPSILON,PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"Complex eigenvalue approximation - not implemented in real scalars"); |
| 491 | } else { | ||
| 492 | 40 | nep->its--; /* do not count this as a full iteration */ | |
| 493 | /* use second eigenpair computed in previous iteration */ | ||
| 494 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(EPSGetConverged(ctx->eps,&nconv)); |
| 495 |
2/4✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
40 | if (nconv>=2 && nconv2>=2) { |
| 496 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(EPSGetEigenpair(ctx->eps,1,&mu,&im,u,NULL)); |
| 497 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(EPSGetEigenpair(ctx->epsts,1,&mu2,&im2,w,NULL)); |
| 498 | 40 | mu = 1.0/mu; | |
| 499 | } else { | ||
| 500 | ✗ | PetscCall(BVSetRandomColumn(nep->V,nep->nconv)); | |
| 501 | ✗ | PetscCall(BVSetRandomColumn(nep->W,nep->nconv)); | |
| 502 | ✗ | PetscCall(BVCopyVec(nep->V,nep->nconv,u)); | |
| 503 | ✗ | PetscCall(BVCopyVec(nep->W,nep->nconv,w)); | |
| 504 | ✗ | mu = lambda-sigma; | |
| 505 | } | ||
| 506 | skip = PETSC_FALSE; | ||
| 507 | } | ||
| 508 | /* correct eigenvalue */ | ||
| 509 | 200 | lambda = lambda - mu; | |
| 510 | } | ||
| 511 | } | ||
| 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.
|
50 | PetscCall(VecDestroy(&u)); |
| 513 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(VecDestroy(&w)); |
| 514 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(VecDestroy(&r)); |
| 515 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(MatDestroy(&mF)); |
| 516 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(MatDestroy(&mJ)); |
| 517 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVDestroy(&X)); |
| 518 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(BVDestroy(&Y)); |
| 519 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
50 | PetscCall(NEPDeflationNEDestroy(defctx)); |
| 520 |
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); |
| 521 | } | ||
| 522 |