| 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 | NEP routines related to the solution process | ||
| 12 | |||
| 13 | References: | ||
| 14 | |||
| 15 | [1] C. Campos and J.E. Roman, "NEP: a module for the parallel solution | ||
| 16 | of nonlinear eigenvalue problems in SLEPc", ACM Trans. Math. Soft. | ||
| 17 | 47(3), 23:1--23:29, 2021. | ||
| 18 | */ | ||
| 19 | |||
| 20 | #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/ | ||
| 21 | #include <slepc/private/bvimpl.h> | ||
| 22 | #include <petscdraw.h> | ||
| 23 | |||
| 24 | static PetscBool cited = PETSC_FALSE; | ||
| 25 | static const char citation[] = | ||
| 26 | "@Article{slepc-nep,\n" | ||
| 27 | " author = \"C. Campos and J. E. Roman\",\n" | ||
| 28 | " title = \"{NEP}: a module for the parallel solution of nonlinear eigenvalue problems in {SLEPc}\",\n" | ||
| 29 | " journal = \"{ACM} Trans. Math. Software\",\n" | ||
| 30 | " volume = \"47\",\n" | ||
| 31 | " number = \"3\",\n" | ||
| 32 | " pages = \"23:1--23:29\",\n" | ||
| 33 | " year = \"2021\",\n" | ||
| 34 | " doi = \"10.1145/3447544\"\n" | ||
| 35 | "}\n"; | ||
| 36 | |||
| 37 | 6844 | PetscErrorCode NEPComputeVectors(NEP nep) | |
| 38 | { | ||
| 39 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
6844 | PetscFunctionBegin; |
| 40 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
6844 | NEPCheckSolved(nep,1); |
| 41 |
8/10✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 8 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
|
6844 | if (nep->state==NEP_STATE_SOLVED) PetscTryTypeMethod(nep,computevectors); |
| 42 | 6844 | nep->state = NEP_STATE_EIGENVECTORS; | |
| 43 |
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.
|
6844 | PetscFunctionReturn(PETSC_SUCCESS); |
| 44 | } | ||
| 45 | |||
| 46 | /*@ | ||
| 47 | NEPSolve - Solves the nonlinear eigenproblem. | ||
| 48 | |||
| 49 | Collective | ||
| 50 | |||
| 51 | Input Parameter: | ||
| 52 | . nep - the nonlinear eigensolver context | ||
| 53 | |||
| 54 | Options Database Keys: | ||
| 55 | + -nep_view - print information about the solver once the solve is complete | ||
| 56 | . -nep_view_pre - print information about the solver before the solve starts | ||
| 57 | . -nep_view_matk - view the split form matrix $A_k$ (replace `k` by an integer from 0 to `nt`-1) | ||
| 58 | . -nep_view_fnk - view the split form function $f_k$ (replace `k` by an integer from 0 to `nt`-1) | ||
| 59 | . -nep_view_vectors - view the computed eigenvectors | ||
| 60 | . -nep_view_values - view the computed eigenvalues | ||
| 61 | . -nep_converged_reason - print reason for convergence/divergence, and number of iterations | ||
| 62 | . -nep_error_absolute - print absolute errors of each eigenpair | ||
| 63 | . -nep_error_relative - print relative errors of each eigenpair | ||
| 64 | - -nep_error_backward - print backward errors of each eigenpair | ||
| 65 | |||
| 66 | Notes: | ||
| 67 | `NEPSolve()` will return without generating an error regardless of whether | ||
| 68 | all requested solutions were computed or not. Call `NEPGetConverged()` to get the | ||
| 69 | actual number of computed solutions, and `NEPGetConvergedReason()` to determine if | ||
| 70 | the solver converged or failed and why. | ||
| 71 | |||
| 72 | All the command-line options listed above admit an optional argument specifying | ||
| 73 | the viewer type and options. For instance, use `-nep_view_vectors binary:myvecs.bin` | ||
| 74 | to save the eigenvectors to a binary file, `-nep_view_values draw` to draw the computed | ||
| 75 | eigenvalues graphically, or `-nep_error_relative :myerr.m:ascii_matlab` to save | ||
| 76 | the errors in a file that can be executed in Matlab. | ||
| 77 | |||
| 78 | Level: beginner | ||
| 79 | |||
| 80 | .seealso: [](ch:nep), `NEPCreate()`, `NEPSetUp()`, `NEPDestroy()`, `NEPSetTolerances()`, `NEPGetConverged()`, `NEPGetConvergedReason()` | ||
| 81 | @*/ | ||
| 82 | 1485 | PetscErrorCode NEPSolve(NEP nep) | |
| 83 | { | ||
| 84 | 1485 | PetscInt i; | |
| 85 | 1485 | char str[16]; | |
| 86 | |||
| 87 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
1485 | PetscFunctionBegin; |
| 88 |
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.
|
1485 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 89 |
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.
|
1485 | if (nep->state>=NEP_STATE_SOLVED) PetscFunctionReturn(PETSC_SUCCESS); |
| 90 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(PetscCitationsRegister(citation,&cited)); |
| 91 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(PetscLogEventBegin(NEP_Solve,nep,0,0,0)); |
| 92 | |||
| 93 | /* call setup */ | ||
| 94 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(NEPSetUp(nep)); |
| 95 | 1485 | nep->nconv = 0; | |
| 96 | 1485 | nep->its = 0; | |
| 97 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
65541 | for (i=0;i<nep->ncv;i++) { |
| 98 | 64056 | nep->eigr[i] = 0.0; | |
| 99 | 64056 | nep->eigi[i] = 0.0; | |
| 100 | 64056 | nep->errest[i] = 0.0; | |
| 101 | 64056 | nep->perm[i] = i; | |
| 102 | } | ||
| 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.
|
1485 | PetscCall(NEPViewFromOptions(nep,NULL,"-nep_view_pre")); |
| 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.
|
1485 | PetscCall(RGViewFromOptions(nep->rg,NULL,"-rg_view")); |
| 105 | |||
| 106 | /* call solver */ | ||
| 107 |
5/10✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✓ Branch 5 taken 8 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
|
1485 | PetscUseTypeMethod(nep,solve); |
| 108 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
1485 | PetscCheck(nep->reason,PetscObjectComm((PetscObject)nep),PETSC_ERR_PLIB,"Internal error, solver returned without setting converged reason"); |
| 109 | 1485 | nep->state = NEP_STATE_SOLVED; | |
| 110 | |||
| 111 | /* Only the first nconv columns contain useful information */ | ||
| 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.
|
1485 | PetscCall(BVSetActiveColumns(nep->V,0,nep->nconv)); |
| 113 |
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.
|
1485 | if (nep->twosided) PetscCall(BVSetActiveColumns(nep->W,0,nep->nconv)); |
| 114 | |||
| 115 |
4/6✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 10 times.
✗ Branch 5 not taken.
|
1485 | if (nep->refine==NEP_REFINE_SIMPLE && nep->rits>0 && nep->nconv>0) { |
| 116 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
136 | PetscCall(NEPComputeVectors(nep)); |
| 117 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
136 | PetscCall(NEPNewtonRefinementSimple(nep,&nep->rits,nep->rtol,nep->nconv)); |
| 118 | 136 | nep->state = NEP_STATE_EIGENVECTORS; | |
| 119 | } | ||
| 120 | |||
| 121 | /* sort eigenvalues according to nep->which parameter */ | ||
| 122 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(SlepcSortEigenvalues(nep->sc,nep->nconv,nep->eigr,nep->eigi,nep->perm)); |
| 123 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(PetscLogEventEnd(NEP_Solve,nep,0,0,0)); |
| 124 | |||
| 125 | /* various viewers */ | ||
| 126 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(NEPViewFromOptions(nep,NULL,"-nep_view")); |
| 127 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(NEPConvergedReasonViewFromOptions(nep)); |
| 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.
|
1485 | PetscCall(NEPErrorViewFromOptions(nep)); |
| 129 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(NEPValuesViewFromOptions(nep)); |
| 130 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(NEPVectorsViewFromOptions(nep)); |
| 131 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1485 | if (nep->fui==NEP_USER_INTERFACE_SPLIT) { |
| 132 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
4016 | for (i=0;i<nep->nt;i++) { |
| 133 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2984 | PetscCall(PetscSNPrintf(str,sizeof(str),"-nep_view_mat%" PetscInt_FMT,i)); |
| 134 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2984 | PetscCall(MatViewFromOptions(nep->A[i],(PetscObject)nep,str)); |
| 135 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2984 | PetscCall(PetscSNPrintf(str,sizeof(str),"-nep_view_fn%" PetscInt_FMT,i)); |
| 136 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2984 | PetscCall(FNViewFromOptions(nep->f[i],(PetscObject)nep,str)); |
| 137 | } | ||
| 138 | } | ||
| 139 | |||
| 140 | /* Remove the initial subspace */ | ||
| 141 | 1485 | nep->nini = 0; | |
| 142 | |||
| 143 | /* Reset resolvent information */ | ||
| 144 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1485 | PetscCall(MatDestroy(&nep->resolvent)); |
| 145 |
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.
|
283 | PetscFunctionReturn(PETSC_SUCCESS); |
| 146 | } | ||
| 147 | |||
| 148 | /*@ | ||
| 149 | NEPProjectOperator - Computes the projection of the nonlinear operator. | ||
| 150 | |||
| 151 | Collective | ||
| 152 | |||
| 153 | Input Parameters: | ||
| 154 | + nep - the nonlinear eigensolver context | ||
| 155 | . j0 - initial index | ||
| 156 | - j1 - final index | ||
| 157 | |||
| 158 | Notes: | ||
| 159 | This is available for split operator only. | ||
| 160 | |||
| 161 | The nonlinear operator $T(\lambda)$ is projected onto $\operatorname{span}(V)$, where $V$ is | ||
| 162 | an orthonormal basis built internally by the solver. The projected | ||
| 163 | operator is equal to $\sum_i V^* A_i V f_i(\lambda)$, so this function | ||
| 164 | computes all matrices $E_i = V^* A_i V$, and stores them in the extra | ||
| 165 | matrices inside `DS`, see `DSMatType`. Only rows/columns in the range [`j0`,`j1`-1] are computed, | ||
| 166 | the previous ones are assumed to be available already. | ||
| 167 | |||
| 168 | Level: developer | ||
| 169 | |||
| 170 | .seealso: [](ch:nep), `NEPSetSplitOperator()` | ||
| 171 | @*/ | ||
| 172 | 10 | PetscErrorCode NEPProjectOperator(NEP nep,PetscInt j0,PetscInt j1) | |
| 173 | { | ||
| 174 | 10 | PetscInt k; | |
| 175 | 10 | Mat G; | |
| 176 | |||
| 177 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBegin; |
| 178 |
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(nep,NEP_CLASSID,1); |
| 179 |
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.
|
10 | PetscValidLogicalCollectiveInt(nep,j0,2); |
| 180 |
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.
|
10 | PetscValidLogicalCollectiveInt(nep,j1,3); |
| 181 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
10 | NEPCheckProblem(nep,1); |
| 182 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
10 | NEPCheckSplit(nep,1); |
| 183 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10 | PetscCall(BVSetActiveColumns(nep->V,j0,j1)); |
| 184 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
40 | for (k=0;k<nep->nt;k++) { |
| 185 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
30 | PetscCall(DSGetMat(nep->ds,DSMatExtra[k],&G)); |
| 186 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
30 | PetscCall(BVMatProject(nep->V,nep->A[k],nep->V,G)); |
| 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.
|
30 | PetscCall(DSRestoreMat(nep->ds,DSMatExtra[k],&G)); |
| 188 | } | ||
| 189 |
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); |
| 190 | } | ||
| 191 | |||
| 192 | /*@ | ||
| 193 | NEPApplyFunction - Applies the nonlinear function $T(\lambda)$ to a given vector. | ||
| 194 | |||
| 195 | Collective | ||
| 196 | |||
| 197 | Input Parameters: | ||
| 198 | + nep - the nonlinear eigensolver context | ||
| 199 | . lambda - scalar argument | ||
| 200 | . x - vector to be multiplied against | ||
| 201 | - v - workspace vector (used only in the case of split form) | ||
| 202 | |||
| 203 | Output Parameters: | ||
| 204 | + y - result vector | ||
| 205 | . A - (optional) Function matrix, for callback interface only | ||
| 206 | - B - (unused) preconditioning matrix | ||
| 207 | |||
| 208 | Note: | ||
| 209 | If the nonlinear operator is represented in split form, the result | ||
| 210 | $y = T(\lambda)x$ is computed without building $T(\lambda)$ explicitly. In | ||
| 211 | that case, parameters `A` and `B` are not used. Otherwise, the matrix | ||
| 212 | $T(\lambda)$ is built and the effect is the same as a call to | ||
| 213 | `NEPComputeFunction()` followed by a `MatMult()`. | ||
| 214 | |||
| 215 | Level: developer | ||
| 216 | |||
| 217 | .seealso: [](ch:nep), `NEPSetSplitOperator()`, `NEPComputeFunction()`, `NEPApplyAdjoint()` | ||
| 218 | @*/ | ||
| 219 | 6928 | PetscErrorCode NEPApplyFunction(NEP nep,PetscScalar lambda,Vec x,Vec v,Vec y,Mat A,Mat B) | |
| 220 | { | ||
| 221 | 6928 | PetscInt i; | |
| 222 | 6928 | PetscScalar alpha; | |
| 223 | |||
| 224 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
6928 | PetscFunctionBegin; |
| 225 |
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.
|
6928 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 226 |
30/68✗ 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 taken 2 times.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 65 taken 2 times.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
|
6928 | PetscValidLogicalCollectiveScalar(nep,lambda,2); |
| 227 |
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.
|
6928 | PetscValidHeaderSpecific(x,VEC_CLASSID,3); |
| 228 |
4/14✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
6928 | if (v) PetscValidHeaderSpecific(v,VEC_CLASSID,4); |
| 229 |
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.
|
6928 | PetscValidHeaderSpecific(y,VEC_CLASSID,5); |
| 230 |
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.
|
6928 | if (A) PetscValidHeaderSpecific(A,MAT_CLASSID,6); |
| 231 |
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.
|
6928 | if (B) PetscValidHeaderSpecific(B,MAT_CLASSID,7); |
| 232 | |||
| 233 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
6928 | if (nep->fui==NEP_USER_INTERFACE_SPLIT) { |
| 234 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5746 | PetscCall(VecSet(y,0.0)); |
| 235 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
22605 | for (i=0;i<nep->nt;i++) { |
| 236 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
16859 | PetscCall(FNEvaluateFunction(nep->f[i],lambda,&alpha)); |
| 237 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
16859 | PetscCall(MatMult(nep->A[i],x,v)); |
| 238 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
16859 | PetscCall(VecAXPY(y,alpha,v)); |
| 239 | } | ||
| 240 | } else { | ||
| 241 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
1182 | if (!A) A = nep->function; |
| 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.
|
1182 | PetscCall(NEPComputeFunction(nep,lambda,A,A)); |
| 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.
|
1182 | PetscCall(MatMult(A,x,y)); |
| 244 | } | ||
| 245 |
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.
|
1361 | PetscFunctionReturn(PETSC_SUCCESS); |
| 246 | } | ||
| 247 | |||
| 248 | /*@ | ||
| 249 | NEPApplyAdjoint - Applies the adjoint nonlinear function $T^*(\lambda)$ to a given vector. | ||
| 250 | |||
| 251 | Collective | ||
| 252 | |||
| 253 | Input Parameters: | ||
| 254 | + nep - the nonlinear eigensolver context | ||
| 255 | . lambda - scalar argument | ||
| 256 | . x - vector to be multiplied against | ||
| 257 | - v - workspace vector (used only in the case of split form) | ||
| 258 | |||
| 259 | Output Parameters: | ||
| 260 | + y - result vector | ||
| 261 | . A - (optional) Function matrix, for callback interface only | ||
| 262 | - B - (unused) preconditioning matrix | ||
| 263 | |||
| 264 | Note: | ||
| 265 | If the nonlinear operator is represented in split form, the result | ||
| 266 | $y = T^*(\lambda)x$ is computed without building $T(\lambda)$ explicitly. In | ||
| 267 | that case, parameters `A` and `B` are not used. Otherwise, the matrix | ||
| 268 | $T(\lambda)$ is built and the effect is the same as a call to | ||
| 269 | `NEPComputeFunction()` followed by a `MatMultHermitianTranspose()`. | ||
| 270 | |||
| 271 | Level: developer | ||
| 272 | |||
| 273 | .seealso: [](ch:nep), `NEPSetSplitOperator()`, `NEPComputeFunction()`, `NEPApplyFunction()` | ||
| 274 | @*/ | ||
| 275 | 210 | PetscErrorCode NEPApplyAdjoint(NEP nep,PetscScalar lambda,Vec x,Vec v,Vec y,Mat A,Mat B) | |
| 276 | { | ||
| 277 | 210 | PetscInt i; | |
| 278 | 210 | PetscScalar alpha; | |
| 279 | 210 | Vec w; | |
| 280 | |||
| 281 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
210 | PetscFunctionBegin; |
| 282 |
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.
|
210 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 283 |
30/68✗ 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 taken 2 times.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 65 taken 2 times.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
|
210 | PetscValidLogicalCollectiveScalar(nep,lambda,2); |
| 284 |
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.
|
210 | PetscValidHeaderSpecific(x,VEC_CLASSID,3); |
| 285 |
4/14✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
210 | if (v) PetscValidHeaderSpecific(v,VEC_CLASSID,4); |
| 286 |
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.
|
210 | PetscValidHeaderSpecific(y,VEC_CLASSID,5); |
| 287 |
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.
|
210 | if (A) PetscValidHeaderSpecific(A,MAT_CLASSID,6); |
| 288 |
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.
|
210 | if (B) PetscValidHeaderSpecific(B,MAT_CLASSID,7); |
| 289 | |||
| 290 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
210 | PetscCall(VecDuplicate(x,&w)); |
| 291 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
210 | PetscCall(VecCopy(x,w)); |
| 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.
|
210 | PetscCall(VecConjugate(w)); |
| 293 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
210 | if (nep->fui==NEP_USER_INTERFACE_SPLIT) { |
| 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.
|
130 | PetscCall(VecSet(y,0.0)); |
| 295 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
520 | for (i=0;i<nep->nt;i++) { |
| 296 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
390 | PetscCall(FNEvaluateFunction(nep->f[i],lambda,&alpha)); |
| 297 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
390 | PetscCall(MatMultHermitianTranspose(nep->A[i],w,v)); |
| 298 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
390 | PetscCall(VecAXPY(y,alpha,v)); |
| 299 | } | ||
| 300 | } else { | ||
| 301 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
80 | if (!A) A = nep->function; |
| 302 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
80 | PetscCall(NEPComputeFunction(nep,lambda,A,A)); |
| 303 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
80 | PetscCall(MatMultHermitianTranspose(A,w,y)); |
| 304 | } | ||
| 305 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
210 | PetscCall(VecDestroy(&w)); |
| 306 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
210 | PetscCall(VecConjugate(y)); |
| 307 |
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.
|
42 | PetscFunctionReturn(PETSC_SUCCESS); |
| 308 | } | ||
| 309 | |||
| 310 | /*@ | ||
| 311 | NEPApplyJacobian - Applies the nonlinear Jacobian $T'(\lambda)$ to a given vector. | ||
| 312 | |||
| 313 | Collective | ||
| 314 | |||
| 315 | Input Parameters: | ||
| 316 | + nep - the nonlinear eigensolver context | ||
| 317 | . lambda - scalar argument | ||
| 318 | . x - vector to be multiplied against | ||
| 319 | - v - workspace vector (used only in the case of split form) | ||
| 320 | |||
| 321 | Output Parameters: | ||
| 322 | + y - result vector | ||
| 323 | - A - (optional) Jacobian matrix, for callback interface only | ||
| 324 | |||
| 325 | Note: | ||
| 326 | If the nonlinear operator is represented in split form, the result | ||
| 327 | $y = T'(\lambda)x$ is computed without building $T'(\lambda)$ explicitly. In | ||
| 328 | that case, parameter `A` is not used. Otherwise, the matrix | ||
| 329 | $T'(\lambda)$ is built and the effect is the same as a call to | ||
| 330 | `NEPComputeJacobian()` followed by a `MatMult()`. | ||
| 331 | |||
| 332 | Level: developer | ||
| 333 | |||
| 334 | .seealso: [](ch:nep), `NEPSetSplitOperator()`, `NEPComputeJacobian()` | ||
| 335 | @*/ | ||
| 336 | 320 | PetscErrorCode NEPApplyJacobian(NEP nep,PetscScalar lambda,Vec x,Vec v,Vec y,Mat A) | |
| 337 | { | ||
| 338 | 320 | PetscInt i; | |
| 339 | 320 | PetscScalar alpha; | |
| 340 | |||
| 341 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
320 | PetscFunctionBegin; |
| 342 |
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.
|
320 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 343 |
30/68✗ 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 taken 2 times.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 65 taken 2 times.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
|
320 | PetscValidLogicalCollectiveScalar(nep,lambda,2); |
| 344 |
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.
|
320 | PetscValidHeaderSpecific(x,VEC_CLASSID,3); |
| 345 |
3/14✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
320 | if (v) PetscValidHeaderSpecific(v,VEC_CLASSID,4); |
| 346 |
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.
|
320 | PetscValidHeaderSpecific(y,VEC_CLASSID,5); |
| 347 |
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.
|
320 | if (A) PetscValidHeaderSpecific(A,MAT_CLASSID,6); |
| 348 | |||
| 349 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
320 | if (nep->fui==NEP_USER_INTERFACE_SPLIT) { |
| 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.
|
320 | PetscCall(VecSet(y,0.0)); |
| 351 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1280 | for (i=0;i<nep->nt;i++) { |
| 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.
|
960 | PetscCall(FNEvaluateDerivative(nep->f[i],lambda,&alpha)); |
| 353 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
960 | PetscCall(MatMult(nep->A[i],x,v)); |
| 354 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
960 | PetscCall(VecAXPY(y,alpha,v)); |
| 355 | } | ||
| 356 | } else { | ||
| 357 | ✗ | if (!A) A = nep->jacobian; | |
| 358 | ✗ | PetscCall(NEPComputeJacobian(nep,lambda,A)); | |
| 359 | ✗ | PetscCall(MatMult(A,x,y)); | |
| 360 | } | ||
| 361 |
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.
|
64 | PetscFunctionReturn(PETSC_SUCCESS); |
| 362 | } | ||
| 363 | |||
| 364 | /*@ | ||
| 365 | NEPGetIterationNumber - Gets the current iteration number. If the | ||
| 366 | call to `NEPSolve()` is complete, then it returns the number of iterations | ||
| 367 | carried out by the solution method. | ||
| 368 | |||
| 369 | Not Collective | ||
| 370 | |||
| 371 | Input Parameter: | ||
| 372 | . nep - the nonlinear eigensolver context | ||
| 373 | |||
| 374 | Output Parameter: | ||
| 375 | . its - number of iterations | ||
| 376 | |||
| 377 | Note: | ||
| 378 | During the $i$-th iteration this call returns $i-1$. If `NEPSolve()` is | ||
| 379 | complete, then parameter `its` contains either the iteration number at | ||
| 380 | which convergence was successfully reached, or failure was detected. | ||
| 381 | Call `NEPGetConvergedReason()` to determine if the solver converged or | ||
| 382 | failed and why. | ||
| 383 | |||
| 384 | Level: intermediate | ||
| 385 | |||
| 386 | .seealso: [](ch:nep), `NEPGetConvergedReason()`, `NEPSetTolerances()` | ||
| 387 | @*/ | ||
| 388 | 99 | PetscErrorCode NEPGetIterationNumber(NEP nep,PetscInt *its) | |
| 389 | { | ||
| 390 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
99 | PetscFunctionBegin; |
| 391 |
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.
|
99 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 392 |
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.
|
99 | PetscAssertPointer(its,2); |
| 393 | 99 | *its = nep->its; | |
| 394 |
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.
|
99 | PetscFunctionReturn(PETSC_SUCCESS); |
| 395 | } | ||
| 396 | |||
| 397 | /*@ | ||
| 398 | NEPGetConverged - Gets the number of converged eigenpairs. | ||
| 399 | |||
| 400 | Not Collective | ||
| 401 | |||
| 402 | Input Parameter: | ||
| 403 | . nep - the nonlinear eigensolver context | ||
| 404 | |||
| 405 | Output Parameter: | ||
| 406 | . nconv - number of converged eigenpairs | ||
| 407 | |||
| 408 | Note: | ||
| 409 | This function should be called after `NEPSolve()` has finished. | ||
| 410 | |||
| 411 | The value `nconv` may be different from the number of requested solutions | ||
| 412 | `nev`, but not larger than `ncv`, see `NEPSetDimensions()`. | ||
| 413 | |||
| 414 | Level: beginner | ||
| 415 | |||
| 416 | .seealso: [](ch:nep), `NEPSetDimensions()`, `NEPSolve()`, `NEPGetEigenpair()` | ||
| 417 | @*/ | ||
| 418 | 121 | PetscErrorCode NEPGetConverged(NEP nep,PetscInt *nconv) | |
| 419 | { | ||
| 420 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
121 | PetscFunctionBegin; |
| 421 |
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.
|
121 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 422 |
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.
|
121 | PetscAssertPointer(nconv,2); |
| 423 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
121 | NEPCheckSolved(nep,1); |
| 424 | 121 | *nconv = nep->nconv; | |
| 425 |
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.
|
121 | PetscFunctionReturn(PETSC_SUCCESS); |
| 426 | } | ||
| 427 | |||
| 428 | /*@ | ||
| 429 | NEPGetConvergedReason - Gets the reason why the `NEPSolve()` iteration was | ||
| 430 | stopped. | ||
| 431 | |||
| 432 | Not Collective | ||
| 433 | |||
| 434 | Input Parameter: | ||
| 435 | . nep - the nonlinear eigensolver context | ||
| 436 | |||
| 437 | Output Parameter: | ||
| 438 | . reason - negative value indicates diverged, positive value converged, see | ||
| 439 | `NEPConvergedReason` for the possible values | ||
| 440 | |||
| 441 | Options Database Key: | ||
| 442 | . -nep_converged_reason - print reason for convergence/divergence, and number of iterations | ||
| 443 | |||
| 444 | Note: | ||
| 445 | If this routine is called before or doing the `NEPSolve()` the value of | ||
| 446 | `NEP_CONVERGED_ITERATING` is returned. | ||
| 447 | |||
| 448 | Level: intermediate | ||
| 449 | |||
| 450 | .seealso: [](ch:nep), `NEPSetTolerances()`, `NEPSolve()`, `NEPConvergedReason` | ||
| 451 | @*/ | ||
| 452 | 16 | PetscErrorCode NEPGetConvergedReason(NEP nep,NEPConvergedReason *reason) | |
| 453 | { | ||
| 454 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
16 | PetscFunctionBegin; |
| 455 |
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.
|
16 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 456 |
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.
|
16 | PetscAssertPointer(reason,2); |
| 457 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
16 | NEPCheckSolved(nep,1); |
| 458 | 16 | *reason = nep->reason; | |
| 459 |
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.
|
16 | PetscFunctionReturn(PETSC_SUCCESS); |
| 460 | } | ||
| 461 | |||
| 462 | /*@ | ||
| 463 | NEPGetEigenpair - Gets the `i`-th solution of the eigenproblem as computed by | ||
| 464 | `NEPSolve()`. The solution consists in both the eigenvalue and the eigenvector. | ||
| 465 | |||
| 466 | Collective | ||
| 467 | |||
| 468 | Input Parameters: | ||
| 469 | + nep - the nonlinear eigensolver context | ||
| 470 | - i - index of the solution | ||
| 471 | |||
| 472 | Output Parameters: | ||
| 473 | + eigr - real part of eigenvalue | ||
| 474 | . eigi - imaginary part of eigenvalue | ||
| 475 | . Vr - real part of eigenvector | ||
| 476 | - Vi - imaginary part of eigenvector | ||
| 477 | |||
| 478 | Notes: | ||
| 479 | It is allowed to pass `NULL` for `Vr` and `Vi`, if the eigenvector is not | ||
| 480 | required. Otherwise, the caller must provide valid `Vec` objects, i.e., | ||
| 481 | they must be created by the calling program with e.g. `MatCreateVecs()`. | ||
| 482 | |||
| 483 | If the eigenvalue is real, then `eigi` and `Vi` are set to zero. If PETSc is | ||
| 484 | configured with complex scalars the eigenvalue is stored | ||
| 485 | directly in `eigr` (`eigi` is set to zero) and the eigenvector in `Vr` (`Vi` | ||
| 486 | is set to zero). In any case, the user can pass `NULL` in any argument that | ||
| 487 | is not required. | ||
| 488 | |||
| 489 | The index `i` should be a value between 0 and `nconv`-1 (see `NEPGetConverged()`). | ||
| 490 | Eigenpairs are indexed according to the ordering criterion established | ||
| 491 | with `NEPSetWhichEigenpairs()`. | ||
| 492 | |||
| 493 | The eigenvector is normalized to have unit norm. | ||
| 494 | |||
| 495 | Level: beginner | ||
| 496 | |||
| 497 | .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConverged()`, `NEPSetWhichEigenpairs()`, `NEPGetLeftEigenvector()` | ||
| 498 | @*/ | ||
| 499 | 6428 | PetscErrorCode NEPGetEigenpair(NEP nep,PetscInt i,PetscScalar *eigr,PetscScalar *eigi,Vec Vr,Vec Vi) | |
| 500 | { | ||
| 501 | 6428 | PetscInt k; | |
| 502 | |||
| 503 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
6428 | PetscFunctionBegin; |
| 504 |
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.
|
6428 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 505 |
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.
|
6428 | PetscValidLogicalCollectiveInt(nep,i,2); |
| 506 |
17/46✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 2 times.
✗ Branch 21 not taken.
✓ 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 taken 2 times.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✓ Branch 35 taken 2 times.
✗ Branch 36 not taken.
✗ 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.
|
6428 | if (Vr) { PetscValidHeaderSpecific(Vr,VEC_CLASSID,5); PetscCheckSameComm(nep,1,Vr,5); } |
| 507 |
17/46✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 2 times.
✗ Branch 21 not taken.
✓ 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 taken 2 times.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✓ Branch 35 taken 2 times.
✗ Branch 36 not taken.
✗ 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.
|
6428 | if (Vi) { PetscValidHeaderSpecific(Vi,VEC_CLASSID,6); PetscCheckSameComm(nep,1,Vi,6); } |
| 508 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
6428 | NEPCheckSolved(nep,1); |
| 509 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
6428 | PetscCheck(i>=0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"The index cannot be negative"); |
| 510 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
6428 | PetscCheck(i<nep->nconv,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"The index can be nconv-1 at most, see NEPGetConverged()"); |
| 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.
|
6428 | PetscCall(NEPComputeVectors(nep)); |
| 513 | 6428 | k = nep->perm[i]; | |
| 514 | |||
| 515 | /* eigenvalue */ | ||
| 516 | #if defined(PETSC_USE_COMPLEX) | ||
| 517 |
1/2✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
|
3661 | if (eigr) *eigr = nep->eigr[k]; |
| 518 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
3661 | if (eigi) *eigi = 0; |
| 519 | #else | ||
| 520 |
1/2✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
|
2767 | if (eigr) *eigr = nep->eigr[k]; |
| 521 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
2767 | if (eigi) *eigi = nep->eigi[k]; |
| 522 | #endif | ||
| 523 | |||
| 524 | /* eigenvector */ | ||
| 525 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
6428 | PetscCall(BV_GetEigenvector(nep->V,k,nep->eigi[k],Vr,Vi)); |
| 526 |
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.
|
1257 | PetscFunctionReturn(PETSC_SUCCESS); |
| 527 | } | ||
| 528 | |||
| 529 | /*@ | ||
| 530 | NEPGetLeftEigenvector - Gets the `i`-th left eigenvector as computed by `NEPSolve()`. | ||
| 531 | |||
| 532 | Collective | ||
| 533 | |||
| 534 | Input Parameters: | ||
| 535 | + nep - the nonlinear eigensolver context | ||
| 536 | - i - index of the solution | ||
| 537 | |||
| 538 | Output Parameters: | ||
| 539 | + Wr - real part of left eigenvector | ||
| 540 | - Wi - imaginary part of left eigenvector | ||
| 541 | |||
| 542 | Notes: | ||
| 543 | The caller must provide valid `Vec` objects, i.e., they must be created | ||
| 544 | by the calling program with e.g. `MatCreateVecs()`. | ||
| 545 | |||
| 546 | If the corresponding eigenvalue is real, then `Wi` is set to zero. If PETSc is | ||
| 547 | configured with complex scalars the eigenvector is stored directly in `Wr` | ||
| 548 | (`Wi` is set to zero). In any case, the user can pass `NULL` in `Wr` or `Wi` if | ||
| 549 | one of them is not required. | ||
| 550 | |||
| 551 | The index `i` should be a value between 0 and `nconv`-1 (see `NEPGetConverged()`). | ||
| 552 | Eigensolutions are indexed according to the ordering criterion established | ||
| 553 | with `NEPSetWhichEigenpairs()`. | ||
| 554 | |||
| 555 | Left eigenvectors are available only if the twosided flag was set, see | ||
| 556 | `NEPSetTwoSided()`. | ||
| 557 | |||
| 558 | Level: intermediate | ||
| 559 | |||
| 560 | .seealso: [](ch:nep), `NEPGetEigenpair()`, `NEPGetConverged()`, `NEPSetWhichEigenpairs()`, `NEPSetTwoSided()` | ||
| 561 | @*/ | ||
| 562 | 210 | PetscErrorCode NEPGetLeftEigenvector(NEP nep,PetscInt i,Vec Wr,Vec Wi) | |
| 563 | { | ||
| 564 | 210 | PetscInt k; | |
| 565 | |||
| 566 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
210 | PetscFunctionBegin; |
| 567 |
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.
|
210 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 568 |
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.
|
210 | PetscValidLogicalCollectiveInt(nep,i,2); |
| 569 |
16/46✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 2 times.
✗ Branch 21 not taken.
✓ 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 taken 2 times.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✓ Branch 35 taken 2 times.
✗ Branch 36 not taken.
✗ 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.
|
210 | if (Wr) { PetscValidHeaderSpecific(Wr,VEC_CLASSID,3); PetscCheckSameComm(nep,1,Wr,3); } |
| 570 |
17/46✓ Branch 0 taken 1 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✓ Branch 14 taken 1 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 1 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 1 times.
✗ Branch 19 not taken.
✓ Branch 20 taken 1 times.
✗ Branch 21 not taken.
✓ Branch 22 taken 1 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 1 times.
✓ Branch 26 taken 1 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 1 times.
✓ Branch 30 taken 1 times.
✗ Branch 31 not taken.
✓ Branch 32 taken 1 times.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✓ Branch 35 taken 1 times.
✗ Branch 36 not taken.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✓ Branch 39 taken 1 times.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✓ Branch 43 taken 1 times.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
|
210 | if (Wi) { PetscValidHeaderSpecific(Wi,VEC_CLASSID,4); PetscCheckSameComm(nep,1,Wi,4); } |
| 571 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
210 | NEPCheckSolved(nep,1); |
| 572 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
210 | PetscCheck(nep->twosided,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONGSTATE,"Must request left vectors with NEPSetTwoSided"); |
| 573 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
210 | PetscCheck(i>=0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"The index cannot be negative"); |
| 574 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
210 | PetscCheck(i<nep->nconv,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"The index can be nconv-1 at most, see NEPGetConverged()"); |
| 575 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
210 | PetscCall(NEPComputeVectors(nep)); |
| 576 | 210 | k = nep->perm[i]; | |
| 577 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
210 | PetscCall(BV_GetEigenvector(nep->W,k,nep->eigi[k],Wr,Wi)); |
| 578 |
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.
|
42 | PetscFunctionReturn(PETSC_SUCCESS); |
| 579 | } | ||
| 580 | |||
| 581 | /*@ | ||
| 582 | NEPGetErrorEstimate - Returns the error estimate associated to the `i`-th | ||
| 583 | computed eigenpair. | ||
| 584 | |||
| 585 | Not Collective | ||
| 586 | |||
| 587 | Input Parameters: | ||
| 588 | + nep - the nonlinear eigensolver context | ||
| 589 | - i - index of eigenpair | ||
| 590 | |||
| 591 | Output Parameter: | ||
| 592 | . errest - the error estimate | ||
| 593 | |||
| 594 | Note: | ||
| 595 | This is the error estimate used internally by the eigensolver. The actual | ||
| 596 | error bound can be computed with `NEPComputeError()`. | ||
| 597 | |||
| 598 | Level: advanced | ||
| 599 | |||
| 600 | .seealso: [](ch:nep), `NEPComputeError()` | ||
| 601 | @*/ | ||
| 602 | 38 | PetscErrorCode NEPGetErrorEstimate(NEP nep,PetscInt i,PetscReal *errest) | |
| 603 | { | ||
| 604 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
38 | PetscFunctionBegin; |
| 605 |
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.
|
38 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 606 |
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.
|
38 | PetscAssertPointer(errest,3); |
| 607 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
38 | NEPCheckSolved(nep,1); |
| 608 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
38 | PetscCheck(i>=0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"The index cannot be negative"); |
| 609 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
38 | PetscCheck(i<nep->nconv,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"The index can be nconv-1 at most, see NEPGetConverged()"); |
| 610 | 38 | *errest = nep->errest[nep->perm[i]]; | |
| 611 |
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.
|
38 | PetscFunctionReturn(PETSC_SUCCESS); |
| 612 | } | ||
| 613 | |||
| 614 | /* | ||
| 615 | NEPComputeResidualNorm_Private - Computes the norm of the residual vector | ||
| 616 | associated with an eigenpair. | ||
| 617 | |||
| 618 | Input Parameters: | ||
| 619 | adj - whether the adjoint T^* must be used instead of T | ||
| 620 | lambda - eigenvalue | ||
| 621 | x - eigenvector | ||
| 622 | w - array of work vectors (two vectors in split form, one vector otherwise) | ||
| 623 | */ | ||
| 624 | 7058 | PetscErrorCode NEPComputeResidualNorm_Private(NEP nep,PetscBool adj,PetscScalar lambda,Vec x,Vec *w,PetscReal *norm) | |
| 625 | { | ||
| 626 | 7058 | Vec y,z=NULL; | |
| 627 | |||
| 628 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
7058 | PetscFunctionBegin; |
| 629 | 7058 | y = w[0]; | |
| 630 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
7058 | if (nep->fui==NEP_USER_INTERFACE_SPLIT) z = w[1]; |
| 631 |
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.
|
7058 | if (adj) PetscCall(NEPApplyAdjoint(nep,lambda,x,z,y,NULL,NULL)); |
| 632 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
6928 | else PetscCall(NEPApplyFunction(nep,lambda,x,z,y,NULL,NULL)); |
| 633 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7058 | PetscCall(VecNorm(y,NORM_2,norm)); |
| 634 |
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.
|
1387 | PetscFunctionReturn(PETSC_SUCCESS); |
| 635 | } | ||
| 636 | |||
| 637 | /*@ | ||
| 638 | NEPComputeError - Computes the error (based on the residual norm) associated | ||
| 639 | with the `i`-th computed eigenpair. | ||
| 640 | |||
| 641 | Collective | ||
| 642 | |||
| 643 | Input Parameters: | ||
| 644 | + nep - the nonlinear eigensolver context | ||
| 645 | . i - the solution index | ||
| 646 | - type - the type of error to compute, see `NEPErrorType` | ||
| 647 | |||
| 648 | Output Parameter: | ||
| 649 | . error - the error | ||
| 650 | |||
| 651 | Notes: | ||
| 652 | The error can be computed in various ways, all of them based on the residual | ||
| 653 | norm computed as $\|T(\lambda)x\|_2$ where $(\lambda,x)$ is the approximate eigenpair. | ||
| 654 | |||
| 655 | If the computation of left eigenvectors was enabled with `NEPSetTwoSided()`, | ||
| 656 | then the error will be computed using the maximum of the value above and | ||
| 657 | the left residual norm $\|y^*T(\lambda)\|_2$, where $y$ is the approximate left | ||
| 658 | eigenvector. | ||
| 659 | |||
| 660 | Level: beginner | ||
| 661 | |||
| 662 | .seealso: [](ch:nep), `NEPErrorType`, `NEPSolve()`, `NEPGetErrorEstimate()`, `NEPSetTwoSided()` | ||
| 663 | @*/ | ||
| 664 | 6091 | PetscErrorCode NEPComputeError(NEP nep,PetscInt i,NEPErrorType type,PetscReal *error) | |
| 665 | { | ||
| 666 | 6091 | Vec xr,xi=NULL; | |
| 667 | 6091 | PetscInt j,nwork,issplit=0; | |
| 668 | 6091 | PetscScalar kr,ki,s; | |
| 669 | 6091 | PetscReal er,z=0.0,errorl,nrm; | |
| 670 | 6091 | PetscBool flg; | |
| 671 | |||
| 672 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
6091 | PetscFunctionBegin; |
| 673 |
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.
|
6091 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 674 |
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.
|
6091 | PetscValidLogicalCollectiveInt(nep,i,2); |
| 675 |
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.
|
6091 | PetscValidLogicalCollectiveEnum(nep,type,3); |
| 676 |
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.
|
6091 | PetscAssertPointer(error,4); |
| 677 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
6091 | NEPCheckSolved(nep,1); |
| 678 | |||
| 679 | /* allocate work vectors */ | ||
| 680 | #if defined(PETSC_USE_COMPLEX) | ||
| 681 | 3488 | nwork = 2; | |
| 682 | #else | ||
| 683 | 2603 | nwork = 3; | |
| 684 | #endif | ||
| 685 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
6091 | if (nep->fui==NEP_USER_INTERFACE_SPLIT) { |
| 686 | 5149 | issplit = 1; | |
| 687 | 5149 | nwork++; /* need an extra work vector for NEPComputeResidualNorm_Private */ | |
| 688 | } | ||
| 689 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
6091 | PetscCall(NEPSetWorkVecs(nep,nwork)); |
| 690 | 6091 | xr = nep->work[issplit+1]; | |
| 691 | #if !defined(PETSC_USE_COMPLEX) | ||
| 692 | 2603 | xi = nep->work[issplit+2]; | |
| 693 | #endif | ||
| 694 | |||
| 695 | /* compute residual norms */ | ||
| 696 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
6091 | PetscCall(NEPGetEigenpair(nep,i,&kr,&ki,xr,xi)); |
| 697 | #if !defined(PETSC_USE_COMPLEX) | ||
| 698 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
2603 | PetscCheck(ki==0.0,PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"Not implemented for complex eigenvalues with real scalars"); |
| 699 | #endif | ||
| 700 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
6091 | PetscCall(NEPComputeResidualNorm_Private(nep,PETSC_FALSE,kr,xr,nep->work,error)); |
| 701 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
6091 | PetscCall(VecNorm(xr,NORM_2,&er)); |
| 702 | |||
| 703 | /* if two-sided, compute left residual norm and take the maximum */ | ||
| 704 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
6091 | if (nep->twosided) { |
| 705 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
130 | PetscCall(NEPGetLeftEigenvector(nep,i,xr,xi)); |
| 706 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
130 | PetscCall(NEPComputeResidualNorm_Private(nep,PETSC_TRUE,kr,xr,nep->work,&errorl)); |
| 707 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
212 | *error = PetscMax(*error,errorl); |
| 708 | } | ||
| 709 | |||
| 710 | /* compute error */ | ||
| 711 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
6091 | switch (type) { |
| 712 | case NEP_ERROR_ABSOLUTE: | ||
| 713 | break; | ||
| 714 | 5425 | case NEP_ERROR_RELATIVE: | |
| 715 | 5425 | *error /= PetscAbsScalar(kr)*er; | |
| 716 | 5425 | break; | |
| 717 | 656 | case NEP_ERROR_BACKWARD: | |
| 718 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
656 | if (nep->fui!=NEP_USER_INTERFACE_SPLIT) { |
| 719 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
270 | PetscCall(NEPComputeFunction(nep,kr,nep->function,nep->function)); |
| 720 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
270 | PetscCall(MatHasOperation(nep->function,MATOP_NORM,&flg)); |
| 721 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
270 | PetscCheck(flg,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"The computation of backward errors requires a matrix norm operation"); |
| 722 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
270 | PetscCall(MatNorm(nep->function,NORM_INFINITY,&nrm)); |
| 723 | 270 | *error /= nrm*er; | |
| 724 | 270 | break; | |
| 725 | } | ||
| 726 | /* initialization of matrix norms */ | ||
| 727 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
386 | if (!nep->nrma[0]) { |
| 728 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
380 | for (j=0;j<nep->nt;j++) { |
| 729 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
266 | PetscCall(MatHasOperation(nep->A[j],MATOP_NORM,&flg)); |
| 730 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
266 | PetscCheck(flg,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"The computation of backward errors requires a matrix norm operation"); |
| 731 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
266 | PetscCall(MatNorm(nep->A[j],NORM_INFINITY,&nep->nrma[j])); |
| 732 | } | ||
| 733 | } | ||
| 734 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1260 | for (j=0;j<nep->nt;j++) { |
| 735 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
874 | PetscCall(FNEvaluateFunction(nep->f[j],kr,&s)); |
| 736 | 874 | z = z + nep->nrma[j]*PetscAbsScalar(s); | |
| 737 | } | ||
| 738 | 386 | *error /= z*er; | |
| 739 | 386 | break; | |
| 740 | ✗ | default: | |
| 741 | ✗ | SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid error type"); | |
| 742 | } | ||
| 743 |
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.
|
1194 | PetscFunctionReturn(PETSC_SUCCESS); |
| 744 | } | ||
| 745 | |||
| 746 | /*@ | ||
| 747 | NEPComputeFunction - Computes the function matrix $T(\lambda)$ that has been | ||
| 748 | set with `NEPSetFunction()`. | ||
| 749 | |||
| 750 | Collective | ||
| 751 | |||
| 752 | Input Parameters: | ||
| 753 | + nep - the nonlinear eigensolver context | ||
| 754 | - lambda - the scalar argument | ||
| 755 | |||
| 756 | Output Parameters: | ||
| 757 | + A - Function matrix | ||
| 758 | - B - optional preconditioning matrix | ||
| 759 | |||
| 760 | Note: | ||
| 761 | `NEPComputeFunction()` is typically used within nonlinear eigensolvers | ||
| 762 | implementations, so most users would not generally call this routine | ||
| 763 | themselves. | ||
| 764 | |||
| 765 | Level: developer | ||
| 766 | |||
| 767 | .seealso: [](ch:nep), `NEPSetFunction()`, `NEPGetFunction()`, `NEPComputeJacobian()` | ||
| 768 | @*/ | ||
| 769 | 29834 | PetscErrorCode NEPComputeFunction(NEP nep,PetscScalar lambda,Mat A,Mat B) | |
| 770 | { | ||
| 771 | 29834 | PetscInt i; | |
| 772 | 29834 | PetscScalar alpha; | |
| 773 | |||
| 774 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
29834 | PetscFunctionBegin; |
| 775 |
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.
|
29834 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 776 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
29834 | NEPCheckProblem(nep,1); |
| 777 |
2/3✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
|
29834 | switch (nep->fui) { |
| 778 | 17995 | case NEP_USER_INTERFACE_CALLBACK: | |
| 779 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
17995 | PetscCheck(nep->computefunction,PetscObjectComm((PetscObject)nep),PETSC_ERR_USER,"Must call NEPSetFunction() first"); |
| 780 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
17995 | PetscCall(PetscLogEventBegin(NEP_FunctionEval,nep,A,B,0)); |
| 781 |
14/26✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 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 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 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
|
17995 | PetscCallBack("NEP user Function function",(*nep->computefunction)(nep,lambda,A,B,nep->functionctx)); |
| 782 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
17995 | PetscCall(PetscLogEventEnd(NEP_FunctionEval,nep,A,B,0)); |
| 783 | break; | ||
| 784 | 11839 | case NEP_USER_INTERFACE_SPLIT: | |
| 785 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
11839 | PetscCall(MatZeroEntries(A)); |
| 786 |
5/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
|
11839 | if (A != B) PetscCall(MatZeroEntries(B)); |
| 787 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
47136 | for (i=0;i<nep->nt;i++) { |
| 788 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35297 | PetscCall(FNEvaluateFunction(nep->f[i],lambda,&alpha)); |
| 789 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35297 | PetscCall(MatAXPY(A,alpha,nep->A[i],nep->mstr)); |
| 790 |
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.
|
35297 | if (A != B) PetscCall(MatAXPY(B,alpha,nep->P[i],nep->mstrp)); |
| 791 | } | ||
| 792 | break; | ||
| 793 | } | ||
| 794 |
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.
|
5680 | PetscFunctionReturn(PETSC_SUCCESS); |
| 795 | } | ||
| 796 | |||
| 797 | /*@ | ||
| 798 | NEPComputeJacobian - Computes the Jacobian matrix $T'(\lambda)$ that has been | ||
| 799 | set with `NEPSetJacobian()`. | ||
| 800 | |||
| 801 | Collective | ||
| 802 | |||
| 803 | Input Parameters: | ||
| 804 | + nep - the nonlinear eigensolver context | ||
| 805 | - lambda - the scalar argument | ||
| 806 | |||
| 807 | Output Parameter: | ||
| 808 | . A - Jacobian matrix | ||
| 809 | |||
| 810 | Notes: | ||
| 811 | Most users should not need to explicitly call this routine, as it | ||
| 812 | is used internally within the nonlinear eigensolvers. | ||
| 813 | |||
| 814 | Level: developer | ||
| 815 | |||
| 816 | .seealso: [](ch:nep), `NEPSetJacobian()`, `NEPGetJacobian()`, `NEPComputeFunction()` | ||
| 817 | @*/ | ||
| 818 | 13846 | PetscErrorCode NEPComputeJacobian(NEP nep,PetscScalar lambda,Mat A) | |
| 819 | { | ||
| 820 | 13846 | PetscInt i; | |
| 821 | 13846 | PetscScalar alpha; | |
| 822 | |||
| 823 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
13846 | PetscFunctionBegin; |
| 824 |
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.
|
13846 | PetscValidHeaderSpecific(nep,NEP_CLASSID,1); |
| 825 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
13846 | NEPCheckProblem(nep,1); |
| 826 |
2/3✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
|
13846 | switch (nep->fui) { |
| 827 | 5120 | case NEP_USER_INTERFACE_CALLBACK: | |
| 828 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
5120 | PetscCheck(nep->computejacobian,PetscObjectComm((PetscObject)nep),PETSC_ERR_USER,"Must call NEPSetJacobian() first"); |
| 829 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5120 | PetscCall(PetscLogEventBegin(NEP_JacobianEval,nep,A,0,0)); |
| 830 |
14/26✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 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 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 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
|
5120 | PetscCallBack("NEP user Jacobian function",(*nep->computejacobian)(nep,lambda,A,nep->jacobianctx)); |
| 831 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5120 | PetscCall(PetscLogEventEnd(NEP_JacobianEval,nep,A,0,0)); |
| 832 | break; | ||
| 833 | 8726 | case NEP_USER_INTERFACE_SPLIT: | |
| 834 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
8726 | PetscCall(MatZeroEntries(A)); |
| 835 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
34736 | for (i=0;i<nep->nt;i++) { |
| 836 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
26010 | PetscCall(FNEvaluateDerivative(nep->f[i],lambda,&alpha)); |
| 837 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
26010 | PetscCall(MatAXPY(A,alpha,nep->A[i],nep->mstr)); |
| 838 | } | ||
| 839 | break; | ||
| 840 | } | ||
| 841 |
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.
|
2684 | PetscFunctionReturn(PETSC_SUCCESS); |
| 842 | } | ||
| 843 |