| 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 | The ST interface routines, callable by users | ||
| 12 | */ | ||
| 13 | |||
| 14 | #include <slepc/private/stimpl.h> /*I "slepcst.h" I*/ | ||
| 15 | |||
| 16 | PetscClassId ST_CLASSID = 0; | ||
| 17 | PetscLogEvent ST_SetUp = 0,ST_ComputeOperator = 0,ST_Apply = 0,ST_ApplyTranspose = 0,ST_ApplyHermitianTranspose = 0,ST_MatSetUp = 0,ST_MatMult = 0,ST_MatMultTranspose = 0,ST_MatSolve = 0,ST_MatSolveTranspose = 0; | ||
| 18 | static PetscBool STPackageInitialized = PETSC_FALSE; | ||
| 19 | |||
| 20 | const char *STMatModes[] = {"COPY","INPLACE","SHELL","STMatMode","ST_MATMODE_",NULL}; | ||
| 21 | |||
| 22 | /*@C | ||
| 23 | STFinalizePackage - This function destroys everything in the SLEPc interface | ||
| 24 | to the `ST` package. It is called from `SlepcFinalize()`. | ||
| 25 | |||
| 26 | Level: developer | ||
| 27 | |||
| 28 | .seealso: [](ch:st), `SlepcFinalize()`, `STInitializePackage()` | ||
| 29 | @*/ | ||
| 30 | 7999 | PetscErrorCode STFinalizePackage(void) | |
| 31 | { | ||
| 32 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
7999 | PetscFunctionBegin; |
| 33 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscFunctionListDestroy(&STList)); |
| 34 | 7999 | STPackageInitialized = PETSC_FALSE; | |
| 35 | 7999 | STRegisterAllCalled = PETSC_FALSE; | |
| 36 |
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.
|
7999 | PetscFunctionReturn(PETSC_SUCCESS); |
| 37 | } | ||
| 38 | |||
| 39 | /*@C | ||
| 40 | STInitializePackage - This function initializes everything in the `ST` package. | ||
| 41 | It is called from `PetscDLLibraryRegister_slepc()` when using dynamic libraries, and | ||
| 42 | on the first call to `STCreate()` when using shared or static libraries. | ||
| 43 | |||
| 44 | Note: | ||
| 45 | This function never needs to be called by SLEPc users. | ||
| 46 | |||
| 47 | Level: developer | ||
| 48 | |||
| 49 | .seealso: [](ch:st), `ST`, `SlepcInitialize()`, `STFinalizePackage()` | ||
| 50 | @*/ | ||
| 51 | 57187 | PetscErrorCode STInitializePackage(void) | |
| 52 | { | ||
| 53 | 57187 | char logList[256]; | |
| 54 | 57187 | PetscBool opt,pkg; | |
| 55 | 57187 | PetscClassId classids[1]; | |
| 56 | |||
| 57 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
57187 | PetscFunctionBegin; |
| 58 |
8/14✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
|
57187 | if (STPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS); |
| 59 | 7999 | STPackageInitialized = PETSC_TRUE; | |
| 60 | /* Register Classes */ | ||
| 61 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscClassIdRegister("Spectral Transform",&ST_CLASSID)); |
| 62 | /* Register Constructors */ | ||
| 63 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(STRegisterAll()); |
| 64 | /* Register Events */ | ||
| 65 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STSetUp",ST_CLASSID,&ST_SetUp)); |
| 66 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STComputeOperatr",ST_CLASSID,&ST_ComputeOperator)); |
| 67 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STApply",ST_CLASSID,&ST_Apply)); |
| 68 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STApplyTranspose",ST_CLASSID,&ST_ApplyTranspose)); |
| 69 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STApplyHermTrans",ST_CLASSID,&ST_ApplyHermitianTranspose)); |
| 70 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STMatSetUp",ST_CLASSID,&ST_MatSetUp)); |
| 71 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STMatMult",ST_CLASSID,&ST_MatMult)); |
| 72 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STMatMultTranspose",ST_CLASSID,&ST_MatMultTranspose)); |
| 73 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STMatSolve",ST_CLASSID,&ST_MatSolve)); |
| 74 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscLogEventRegister("STMatSolveTranspose",ST_CLASSID,&ST_MatSolveTranspose)); |
| 75 | /* Process Info */ | ||
| 76 | 7999 | classids[0] = ST_CLASSID; | |
| 77 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscInfoProcessClass("st",1,&classids[0])); |
| 78 | /* Process summary exclusions */ | ||
| 79 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt)); |
| 80 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
7999 | if (opt) { |
| 81 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
84 | PetscCall(PetscStrInList("st",logList,',',&pkg)); |
| 82 |
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.
|
84 | if (pkg) PetscCall(PetscLogEventDeactivateClass(ST_CLASSID)); |
| 83 | } | ||
| 84 | /* Register package finalizer */ | ||
| 85 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
7999 | PetscCall(PetscRegisterFinalize(STFinalizePackage)); |
| 86 |
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.
|
1571 | PetscFunctionReturn(PETSC_SUCCESS); |
| 87 | } | ||
| 88 | |||
| 89 | /*@ | ||
| 90 | STReset - Resets the `ST` context to the initial state (prior to setup) | ||
| 91 | and destroys any allocated `Vec`s and `Mat`s. | ||
| 92 | |||
| 93 | Collective | ||
| 94 | |||
| 95 | Input Parameter: | ||
| 96 | . st - the spectral transformation context | ||
| 97 | |||
| 98 | Level: advanced | ||
| 99 | |||
| 100 | .seealso: [](ch:st), `STDestroy()` | ||
| 101 | @*/ | ||
| 102 | 20368 | PetscErrorCode STReset(ST st) | |
| 103 | { | ||
| 104 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
20368 | PetscFunctionBegin; |
| 105 |
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.
|
20368 | if (st) PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 106 |
2/14✓ Branch 0 taken 8 times.
✓ Branch 1 taken 1 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.
|
18357 | if (!st) PetscFunctionReturn(PETSC_SUCCESS); |
| 107 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
20368 | STCheckNotSeized(st,1); |
| 108 |
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.
|
20368 | PetscTryTypeMethod(st,reset); |
| 109 |
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.
|
20368 | if (st->ksp) PetscCall(KSPReset(st->ksp)); |
| 110 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20368 | PetscCall(MatDestroyMatrices(PetscMax(2,st->nmat),&st->T)); |
| 111 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20368 | PetscCall(MatDestroyMatrices(PetscMax(2,st->nmat),&st->A)); |
| 112 | 20368 | st->nmat = 0; | |
| 113 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
20368 | PetscCall(PetscFree(st->Astate)); |
| 114 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20368 | PetscCall(MatDestroy(&st->Op)); |
| 115 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20368 | PetscCall(MatDestroy(&st->P)); |
| 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.
|
20368 | PetscCall(MatDestroy(&st->Pmat)); |
| 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.
|
20368 | PetscCall(MatDestroyMatrices(st->nsplit,&st->Psplit)); |
| 118 | 20368 | st->nsplit = 0; | |
| 119 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20368 | PetscCall(VecDestroyVecs(st->nwork,&st->work)); |
| 120 | 20368 | st->nwork = 0; | |
| 121 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20368 | PetscCall(VecDestroy(&st->wb)); |
| 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.
|
20368 | PetscCall(VecDestroy(&st->wht)); |
| 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.
|
20368 | PetscCall(VecDestroy(&st->D)); |
| 124 | 20368 | st->state = ST_STATE_INITIAL; | |
| 125 | 20368 | st->opready = PETSC_FALSE; | |
| 126 |
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.
|
20368 | PetscFunctionReturn(PETSC_SUCCESS); |
| 127 | } | ||
| 128 | |||
| 129 | /*@ | ||
| 130 | STDestroy - Destroys ST context that was created with `STCreate()`. | ||
| 131 | |||
| 132 | Collective | ||
| 133 | |||
| 134 | Input Parameter: | ||
| 135 | . st - the spectral transformation context | ||
| 136 | |||
| 137 | Level: beginner | ||
| 138 | |||
| 139 | .seealso: [](ch:st), `STCreate()`, `STSetUp()` | ||
| 140 | @*/ | ||
| 141 | 9393 | PetscErrorCode STDestroy(ST *st) | |
| 142 | { | ||
| 143 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
9393 | PetscFunctionBegin; |
| 144 |
8/14✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
|
9393 | if (!*st) PetscFunctionReturn(PETSC_SUCCESS); |
| 145 |
2/12✗ 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 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
9189 | PetscValidHeaderSpecific(*st,ST_CLASSID,1); |
| 146 |
8/14✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
|
9189 | if (--((PetscObject)*st)->refct > 0) { *st = NULL; PetscFunctionReturn(PETSC_SUCCESS); } |
| 147 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9148 | PetscCall(STReset(*st)); |
| 148 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
9148 | PetscTryTypeMethod(*st,destroy); |
| 149 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9148 | PetscCall(KSPDestroy(&(*st)->ksp)); |
| 150 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9148 | PetscCall(PetscHeaderDestroy(st)); |
| 151 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
1765 | PetscFunctionReturn(PETSC_SUCCESS); |
| 152 | } | ||
| 153 | |||
| 154 | /*@ | ||
| 155 | STCreate - Creates a spectral transformation context. | ||
| 156 | |||
| 157 | Collective | ||
| 158 | |||
| 159 | Input Parameter: | ||
| 160 | . comm - MPI communicator | ||
| 161 | |||
| 162 | Output Parameter: | ||
| 163 | . newst - location to put the spectral transformation context | ||
| 164 | |||
| 165 | Level: beginner | ||
| 166 | |||
| 167 | .seealso: [](ch:st), `STSetUp()`, `STApply()`, `STDestroy()`, `ST` | ||
| 168 | @*/ | ||
| 169 | 9148 | PetscErrorCode STCreate(MPI_Comm comm,ST *newst) | |
| 170 | { | ||
| 171 | 9148 | ST st; | |
| 172 | |||
| 173 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
9148 | PetscFunctionBegin; |
| 174 |
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.
|
9148 | PetscAssertPointer(newst,2); |
| 175 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9148 | PetscCall(STInitializePackage()); |
| 176 |
7/12✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 10 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 8 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
9148 | PetscCall(SlepcHeaderCreate(st,ST_CLASSID,"ST","Spectral Transformation","ST",comm,STDestroy,STView)); |
| 177 | |||
| 178 | 9148 | st->A = NULL; | |
| 179 | 9148 | st->nmat = 0; | |
| 180 | 9148 | st->sigma = 0.0; | |
| 181 | 9148 | st->defsigma = 0.0; | |
| 182 | 9148 | st->matmode = ST_MATMODE_COPY; | |
| 183 | 9148 | st->str = UNKNOWN_NONZERO_PATTERN; | |
| 184 | 9148 | st->transform = PETSC_FALSE; | |
| 185 | 9148 | st->structured = PETSC_FALSE; | |
| 186 | 9148 | st->D = NULL; | |
| 187 | 9148 | st->Pmat = NULL; | |
| 188 | 9148 | st->Pmat_set = PETSC_FALSE; | |
| 189 | 9148 | st->Psplit = NULL; | |
| 190 | 9148 | st->nsplit = 0; | |
| 191 | 9148 | st->strp = UNKNOWN_NONZERO_PATTERN; | |
| 192 | |||
| 193 | 9148 | st->ksp = NULL; | |
| 194 | 9148 | st->usesksp = PETSC_FALSE; | |
| 195 | 9148 | st->nwork = 0; | |
| 196 | 9148 | st->work = NULL; | |
| 197 | 9148 | st->wb = NULL; | |
| 198 | 9148 | st->wht = NULL; | |
| 199 | 9148 | st->state = ST_STATE_INITIAL; | |
| 200 | 9148 | st->Astate = NULL; | |
| 201 | 9148 | st->T = NULL; | |
| 202 | 9148 | st->Op = NULL; | |
| 203 | 9148 | st->opseized = PETSC_FALSE; | |
| 204 | 9148 | st->opready = PETSC_FALSE; | |
| 205 | 9148 | st->P = NULL; | |
| 206 | 9148 | st->M = NULL; | |
| 207 | 9148 | st->sigma_set = PETSC_FALSE; | |
| 208 | 9148 | st->asymm = PETSC_FALSE; | |
| 209 | 9148 | st->aherm = PETSC_FALSE; | |
| 210 | 9148 | st->data = NULL; | |
| 211 | |||
| 212 | 9148 | *newst = st; | |
| 213 |
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.
|
9148 | PetscFunctionReturn(PETSC_SUCCESS); |
| 214 | } | ||
| 215 | |||
| 216 | /* | ||
| 217 | Checks whether the ST matrices are all symmetric or hermitian. | ||
| 218 | */ | ||
| 219 | 9838 | static inline PetscErrorCode STMatIsSymmetricKnown(ST st,PetscBool *symm,PetscBool *herm) | |
| 220 | { | ||
| 221 | 9838 | PetscInt i; | |
| 222 | 9838 | PetscBool sbaij=PETSC_FALSE,set,flg=PETSC_FALSE; | |
| 223 | |||
| 224 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
9838 | PetscFunctionBegin; |
| 225 | /* check if problem matrices are all sbaij */ | ||
| 226 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
9886 | for (i=0;i<st->nmat;i++) { |
| 227 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9862 | PetscCall(PetscObjectTypeCompareAny((PetscObject)st->A[i],&sbaij,MATSEQSBAIJ,MATMPISBAIJ,"")); |
| 228 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
9862 | if (!sbaij) break; |
| 229 | } | ||
| 230 | /* check if user has set the symmetric flag */ | ||
| 231 | 9838 | *symm = PETSC_TRUE; | |
| 232 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10496 | for (i=0;i<st->nmat;i++) { |
| 233 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10168 | PetscCall(MatIsSymmetricKnown(st->A[i],&set,&flg)); |
| 234 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
10168 | if (!set || !flg) { *symm = PETSC_FALSE; break; } |
| 235 | } | ||
| 236 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
9838 | if (sbaij) *symm = PETSC_TRUE; |
| 237 | #if defined(PETSC_USE_COMPLEX) | ||
| 238 | /* check if user has set the hermitian flag */ | ||
| 239 | 4850 | *herm = PETSC_TRUE; | |
| 240 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
5182 | for (i=0;i<st->nmat;i++) { |
| 241 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
4994 | PetscCall(MatIsHermitianKnown(st->A[i],&set,&flg)); |
| 242 |
4/4✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 5 times.
✓ Branch 3 taken 5 times.
|
4994 | if (!set || !flg) { *herm = PETSC_FALSE; break; } |
| 243 | } | ||
| 244 | #else | ||
| 245 | 4988 | *herm = *symm; | |
| 246 | #endif | ||
| 247 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
5938 | PetscFunctionReturn(PETSC_SUCCESS); |
| 248 | } | ||
| 249 | |||
| 250 | /*@ | ||
| 251 | STSetMatrices - Sets the matrices associated with the eigenvalue problem. | ||
| 252 | |||
| 253 | Collective | ||
| 254 | |||
| 255 | Input Parameters: | ||
| 256 | + st - the spectral transformation context | ||
| 257 | . n - number of matrices in array `A` | ||
| 258 | - A - the array of matrices associated with the eigensystem | ||
| 259 | |||
| 260 | Notes: | ||
| 261 | It must be called before `STSetUp()`. If it is called again after `STSetUp()` then | ||
| 262 | the `ST` object is reset. | ||
| 263 | |||
| 264 | In standard eigenproblems only one matrix is passed, while in generalized | ||
| 265 | problems two matrices are provided. The number of matrices is larger in | ||
| 266 | polynomial eigenproblems. | ||
| 267 | |||
| 268 | In normal usage, matrices are provided via the corresponding `EPS` of `PEP` | ||
| 269 | interface function. | ||
| 270 | |||
| 271 | Level: intermediate | ||
| 272 | |||
| 273 | .seealso: [](ch:st), `STGetMatrix()`, `STGetNumMatrices()`, `STSetUp()`, `STReset()` | ||
| 274 | @*/ | ||
| 275 | 10425 | PetscErrorCode STSetMatrices(ST st,PetscInt n,Mat A[]) | |
| 276 | { | ||
| 277 | 10425 | PetscInt i; | |
| 278 | 10425 | PetscBool same=PETSC_TRUE; | |
| 279 | |||
| 280 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10425 | PetscFunctionBegin; |
| 281 |
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.
|
10425 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 282 |
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.
|
10425 | PetscValidLogicalCollectiveInt(st,n,2); |
| 283 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
10425 | PetscCheck(n>0,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more matrices, you have %" PetscInt_FMT,n); |
| 284 |
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.
|
10425 | PetscAssertPointer(A,3); |
| 285 |
13/32✓ 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 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.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
10425 | PetscCheckSameComm(st,1,*A,3); |
| 286 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
10425 | STCheckNotSeized(st,1); |
| 287 |
3/6✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
10425 | PetscCheck(!st->nsplit || st->nsplit==n,PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"The number of matrices must be the same as in STSetSplitPreconditioner()"); |
| 288 | |||
| 289 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10425 | if (st->state) { |
| 290 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
946 | if (n!=st->nmat) same = PETSC_FALSE; |
| 291 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2046 | for (i=0;same&&i<n;i++) { |
| 292 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1100 | if (A[i]!=st->A[i]) same = PETSC_FALSE; |
| 293 | } | ||
| 294 |
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.
|
946 | if (!same) PetscCall(STReset(st)); |
| 295 | } else same = PETSC_FALSE; | ||
| 296 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
|
160 | if (!same) { |
| 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.
|
9838 | PetscCall(MatDestroyMatrices(PetscMax(2,st->nmat),&st->A)); |
| 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.
|
9838 | PetscCall(PetscCalloc1(PetscMax(2,n),&st->A)); |
| 299 |
5/8✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
9838 | PetscCall(PetscFree(st->Astate)); |
| 300 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9838 | PetscCall(PetscMalloc1(PetscMax(2,n),&st->Astate)); |
| 301 | } | ||
| 302 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
28811 | for (i=0;i<n;i++) { |
| 303 |
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.
|
18386 | PetscValidHeaderSpecific(A[i],MAT_CLASSID,3); |
| 304 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
18386 | PetscCall(PetscObjectReference((PetscObject)A[i])); |
| 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.
|
18386 | PetscCall(MatDestroy(&st->A[i])); |
| 306 | 18386 | st->A[i] = A[i]; | |
| 307 | 18386 | st->Astate[i] = ((PetscObject)A[i])->state; | |
| 308 | } | ||
| 309 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10425 | if (n==1) { |
| 310 | 6028 | st->A[1] = NULL; | |
| 311 | 6028 | st->Astate[1] = 0; | |
| 312 | } | ||
| 313 | 10425 | st->nmat = n; | |
| 314 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
10425 | if (same) st->state = ST_STATE_UPDATED; |
| 315 | 9838 | else st->state = ST_STATE_INITIAL; | |
| 316 |
3/6✓ Branch 0 taken 5 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 5 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
10425 | PetscCheck(!same || !st->Psplit,PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Support for changing the matrices while using a split preconditioner is not implemented yet"); |
| 317 | 10425 | st->opready = PETSC_FALSE; | |
| 318 |
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.
|
10425 | if (!same) PetscCall(STMatIsSymmetricKnown(st,&st->asymm,&st->aherm)); |
| 319 |
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.
|
2017 | PetscFunctionReturn(PETSC_SUCCESS); |
| 320 | } | ||
| 321 | |||
| 322 | /*@ | ||
| 323 | STGetMatrix - Gets the matrices associated with the original eigensystem. | ||
| 324 | |||
| 325 | Collective | ||
| 326 | |||
| 327 | Input Parameters: | ||
| 328 | + st - the spectral transformation context | ||
| 329 | - k - the index of the requested matrix (starting in 0) | ||
| 330 | |||
| 331 | Output Parameter: | ||
| 332 | . A - the requested matrix | ||
| 333 | |||
| 334 | Level: intermediate | ||
| 335 | |||
| 336 | .seealso: [](ch:st), `STSetMatrices()`, `STGetNumMatrices()` | ||
| 337 | @*/ | ||
| 338 | 111558 | PetscErrorCode STGetMatrix(ST st,PetscInt k,Mat *A) | |
| 339 | { | ||
| 340 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
111558 | PetscFunctionBegin; |
| 341 |
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.
|
111558 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 342 |
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.
|
111558 | PetscValidLogicalCollectiveInt(st,k,2); |
| 343 |
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.
|
111558 | PetscAssertPointer(A,3); |
| 344 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
111558 | STCheckMatrices(st,1); |
| 345 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
111558 | PetscCheck(k>=0 && k<st->nmat,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %" PetscInt_FMT,st->nmat-1); |
| 346 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
111558 | PetscCheck(((PetscObject)st->A[k])->state==st->Astate[k],PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Cannot retrieve original matrices (have been modified)"); |
| 347 | 111558 | *A = st->A[k]; | |
| 348 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
111558 | PetscFunctionReturn(PETSC_SUCCESS); |
| 349 | } | ||
| 350 | |||
| 351 | /*@ | ||
| 352 | STGetMatrixTransformed - Gets the matrices associated with the transformed eigensystem. | ||
| 353 | |||
| 354 | Not Collective | ||
| 355 | |||
| 356 | Input Parameters: | ||
| 357 | + st - the spectral transformation context | ||
| 358 | - k - the index of the requested matrix (starting in 0) | ||
| 359 | |||
| 360 | Output Parameter: | ||
| 361 | . T - the requested matrix | ||
| 362 | |||
| 363 | Level: developer | ||
| 364 | |||
| 365 | .seealso: [](ch:st), `STGetMatrix()`, `STGetNumMatrices()` | ||
| 366 | @*/ | ||
| 367 | 2379 | PetscErrorCode STGetMatrixTransformed(ST st,PetscInt k,Mat *T) | |
| 368 | { | ||
| 369 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2379 | PetscFunctionBegin; |
| 370 |
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.
|
2379 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 371 |
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.
|
2379 | PetscValidLogicalCollectiveInt(st,k,2); |
| 372 |
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.
|
2379 | PetscAssertPointer(T,3); |
| 373 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
2379 | STCheckMatrices(st,1); |
| 374 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
2379 | PetscCheck(k>=0 && k<st->nmat,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %" PetscInt_FMT,st->nmat-1); |
| 375 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
2379 | PetscCheck(st->T,PetscObjectComm((PetscObject)st),PETSC_ERR_POINTER,"There are no transformed matrices"); |
| 376 | 2379 | *T = st->T[k]; | |
| 377 |
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.
|
2379 | PetscFunctionReturn(PETSC_SUCCESS); |
| 378 | } | ||
| 379 | |||
| 380 | /*@ | ||
| 381 | STGetNumMatrices - Returns the number of matrices stored in the `ST`. | ||
| 382 | |||
| 383 | Not Collective | ||
| 384 | |||
| 385 | Input Parameter: | ||
| 386 | . st - the spectral transformation context | ||
| 387 | |||
| 388 | Output Parameter: | ||
| 389 | . n - the number of matrices passed in `STSetMatrices()` | ||
| 390 | |||
| 391 | Level: intermediate | ||
| 392 | |||
| 393 | .seealso: [](ch:st), `STSetMatrices()` | ||
| 394 | @*/ | ||
| 395 | 72608 | PetscErrorCode STGetNumMatrices(ST st,PetscInt *n) | |
| 396 | { | ||
| 397 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
72608 | PetscFunctionBegin; |
| 398 |
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.
|
72608 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 399 |
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.
|
72608 | PetscAssertPointer(n,2); |
| 400 | 72608 | *n = st->nmat; | |
| 401 |
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.
|
72608 | PetscFunctionReturn(PETSC_SUCCESS); |
| 402 | } | ||
| 403 | |||
| 404 | /*@ | ||
| 405 | STResetMatrixState - Resets the stored state of the matrices in the `ST`. | ||
| 406 | |||
| 407 | Logically Collective | ||
| 408 | |||
| 409 | Input Parameter: | ||
| 410 | . st - the spectral transformation context | ||
| 411 | |||
| 412 | Note: | ||
| 413 | This is useful in solvers where the user matrices are modified during | ||
| 414 | the computation, as in nonlinear inverse iteration. The effect is that | ||
| 415 | `STGetMatrix()` will retrieve the modified matrices as if they were | ||
| 416 | the matrices originally provided by the user. | ||
| 417 | |||
| 418 | Level: developer | ||
| 419 | |||
| 420 | .seealso: [](ch:st), `STGetMatrix()`, `EPSPowerSetNonlinear()` | ||
| 421 | @*/ | ||
| 422 | 14435 | PetscErrorCode STResetMatrixState(ST st) | |
| 423 | { | ||
| 424 | 14435 | PetscInt i; | |
| 425 | |||
| 426 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
14435 | PetscFunctionBegin; |
| 427 |
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.
|
14435 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 428 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
43305 | for (i=0;i<st->nmat;i++) st->Astate[i] = ((PetscObject)st->A[i])->state; |
| 429 |
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.
|
14435 | PetscFunctionReturn(PETSC_SUCCESS); |
| 430 | } | ||
| 431 | |||
| 432 | /*@ | ||
| 433 | STSetPreconditionerMat - Sets the matrix to be used to build the preconditioner. | ||
| 434 | |||
| 435 | Collective | ||
| 436 | |||
| 437 | Input Parameters: | ||
| 438 | + st - the spectral transformation context | ||
| 439 | - mat - the matrix that will be used in constructing the preconditioner | ||
| 440 | |||
| 441 | Notes: | ||
| 442 | This matrix will be passed to the internal `KSP` object (via the last argument | ||
| 443 | of `KSPSetOperators()`) as the matrix to be used when constructing the preconditioner. | ||
| 444 | If no matrix is set or `mat` is set to `NULL`, then $A-\sigma B$ will be used | ||
| 445 | to build the preconditioner, being $\sigma$ the value set by `STSetShift()`. | ||
| 446 | |||
| 447 | More precisely, this is relevant for spectral transformations that represent | ||
| 448 | a rational matrix function, and use a `KSP` object for the denominator, called | ||
| 449 | $K$ in the description of `STGetOperator()`. It includes also the `STPRECOND` case. | ||
| 450 | If the user has a good approximation to matrix $K$ that can be used to build a | ||
| 451 | cheap preconditioner, it can be passed with this function. Note that it affects | ||
| 452 | only the `Pmat` argument of `KSPSetOperators()`, not the `Amat` argument. | ||
| 453 | |||
| 454 | If a preconditioner matrix is set, the default is to use an iterative `KSP` | ||
| 455 | rather than a direct method. | ||
| 456 | |||
| 457 | An alternative to pass an approximation of $A-\sigma B$ with this function is | ||
| 458 | to provide approximations of $A$ and $B$ via `STSetSplitPreconditioner()`. The | ||
| 459 | difference is that when $\sigma$ changes the preconditioner is recomputed. | ||
| 460 | |||
| 461 | Use `NULL` to remove a previously set matrix. | ||
| 462 | |||
| 463 | Level: advanced | ||
| 464 | |||
| 465 | .seealso: [](ch:st), `STGetPreconditionerMat()`, `STSetShift()`, `STGetOperator()`, `STSetSplitPreconditioner()` | ||
| 466 | @*/ | ||
| 467 | 210 | PetscErrorCode STSetPreconditionerMat(ST st,Mat mat) | |
| 468 | { | ||
| 469 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
210 | PetscFunctionBegin; |
| 470 |
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(st,ST_CLASSID,1); |
| 471 |
2/2✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
|
210 | if (mat) { |
| 472 |
2/12✗ 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 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
208 | PetscValidHeaderSpecific(mat,MAT_CLASSID,2); |
| 473 |
13/32✓ 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 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.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
208 | PetscCheckSameComm(st,1,mat,2); |
| 474 | } | ||
| 475 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
210 | STCheckNotSeized(st,1); |
| 476 |
3/6✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
210 | PetscCheck(!mat || !st->Psplit,PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Cannot call both STSetPreconditionerMat and STSetSplitPreconditioner"); |
| 477 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 5 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.
|
205 | if (mat) PetscCall(PetscObjectReference((PetscObject)mat)); |
| 478 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
210 | PetscCall(MatDestroy(&st->Pmat)); |
| 479 | 210 | st->Pmat = mat; | |
| 480 | 210 | st->Pmat_set = mat? PETSC_TRUE: PETSC_FALSE; | |
| 481 | 210 | st->state = ST_STATE_INITIAL; | |
| 482 | 210 | st->opready = PETSC_FALSE; | |
| 483 |
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.
|
210 | PetscFunctionReturn(PETSC_SUCCESS); |
| 484 | } | ||
| 485 | |||
| 486 | /*@ | ||
| 487 | STGetPreconditionerMat - Returns the matrix previously set by `STSetPreconditionerMat()`. | ||
| 488 | |||
| 489 | Not Collective | ||
| 490 | |||
| 491 | Input Parameter: | ||
| 492 | . st - the spectral transformation context | ||
| 493 | |||
| 494 | Output Parameter: | ||
| 495 | . mat - the matrix that will be used in constructing the preconditioner or | ||
| 496 | `NULL` if no matrix was set by `STSetPreconditionerMat()`. | ||
| 497 | |||
| 498 | Level: advanced | ||
| 499 | |||
| 500 | .seealso: [](ch:st), `STSetPreconditionerMat()` | ||
| 501 | @*/ | ||
| 502 | 460 | PetscErrorCode STGetPreconditionerMat(ST st,Mat *mat) | |
| 503 | { | ||
| 504 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
460 | PetscFunctionBegin; |
| 505 |
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.
|
460 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 506 |
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.
|
460 | PetscAssertPointer(mat,2); |
| 507 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
460 | *mat = st->Pmat_set? st->Pmat: NULL; |
| 508 |
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.
|
460 | PetscFunctionReturn(PETSC_SUCCESS); |
| 509 | } | ||
| 510 | |||
| 511 | /*@ | ||
| 512 | STSetSplitPreconditioner - Sets the matrices from which to build the preconditioner | ||
| 513 | in split form. | ||
| 514 | |||
| 515 | Collective | ||
| 516 | |||
| 517 | Input Parameters: | ||
| 518 | + st - the spectral transformation context | ||
| 519 | . n - number of matrices | ||
| 520 | . Psplit - array of matrices | ||
| 521 | - strp - structure flag for `Psplit` matrices | ||
| 522 | |||
| 523 | Notes: | ||
| 524 | The number of matrices passed here must be the same as in `STSetMatrices()`. | ||
| 525 | |||
| 526 | For linear eigenproblems, the preconditioner matrix is computed as | ||
| 527 | $P(\sigma) = A_0-\sigma B_0$, where $A_0$ and $B_0$ are approximations of $A$ and $B$ | ||
| 528 | (the eigenproblem matrices) provided via the `Psplit` array in this function. | ||
| 529 | Compared to `STSetPreconditionerMat()`, this function allows setting a preconditioner | ||
| 530 | in a way that is independent of the shift $\sigma$. Whenever the value of $\sigma$ | ||
| 531 | changes the preconditioner is recomputed. | ||
| 532 | |||
| 533 | Similarly, for polynomial eigenproblems the matrix for the preconditioner | ||
| 534 | is expressed as $P(\sigma) = \sum_i P_i \phi_i(\sigma)$, for $i=1,\dots,n$, where | ||
| 535 | $P_i$ are given in `Psplit` and the $\phi_i$'s are the polynomial basis functions. | ||
| 536 | |||
| 537 | The structure flag provides information about the relative nonzero pattern of the | ||
| 538 | `Psplit` matrices, in the same way as in `STSetMatStructure()`. | ||
| 539 | |||
| 540 | Use `n=0` to reset a previously set split preconditioner. | ||
| 541 | |||
| 542 | Level: advanced | ||
| 543 | |||
| 544 | .seealso: [](ch:st), `STGetSplitPreconditionerTerm()`, `STGetSplitPreconditionerInfo()`, `STSetPreconditionerMat()`, `STSetMatrices()`, `STSetMatStructure()` | ||
| 545 | @*/ | ||
| 546 | 215 | PetscErrorCode STSetSplitPreconditioner(ST st,PetscInt n,Mat Psplit[],MatStructure strp) | |
| 547 | { | ||
| 548 | 215 | PetscInt i,N=0,M,M0=0,mloc,nloc,mloc0=0; | |
| 549 | |||
| 550 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
215 | PetscFunctionBegin; |
| 551 |
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.
|
215 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 552 |
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.
|
215 | PetscValidLogicalCollectiveInt(st,n,2); |
| 553 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
215 | PetscCheck(n>=0,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Negative value of n = %" PetscInt_FMT,n); |
| 554 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
215 | PetscCheck(!n || !st->Pmat_set,PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Cannot call both STSetPreconditionerMat and STSetSplitPreconditioner"); |
| 555 |
6/8✓ Branch 0 taken 9 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 5 times.
✓ Branch 4 taken 3 times.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
215 | PetscCheck(!n || !st->nmat || st->nmat==n,PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"The number of matrices must be the same as in STSetMatrices()"); |
| 556 |
3/10✓ 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.
|
43 | if (n) PetscAssertPointer(Psplit,3); |
| 557 |
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.
|
43 | PetscValidLogicalCollectiveEnum(st,strp,4); |
| 558 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
43 | STCheckNotSeized(st,1); |
| 559 | |||
| 560 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
700 | for (i=0;i<n;i++) { |
| 561 |
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.
|
485 | PetscValidHeaderSpecific(Psplit[i],MAT_CLASSID,3); |
| 562 |
13/32✓ 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 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.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
485 | PetscCheckSameComm(st,1,Psplit[i],3); |
| 563 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
485 | PetscCall(MatGetSize(Psplit[i],&M,&N)); |
| 564 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
485 | PetscCall(MatGetLocalSize(Psplit[i],&mloc,&nloc)); |
| 565 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
485 | PetscCheck(M==N,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_WRONG,"Psplit[%" PetscInt_FMT "] is a non-square matrix (%" PetscInt_FMT " rows, %" PetscInt_FMT " cols)",i,M,N); |
| 566 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
485 | PetscCheck(mloc==nloc,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_WRONG,"Psplit[%" PetscInt_FMT "] does not have equal row and column local sizes (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,mloc,nloc); |
| 567 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
485 | if (!i) { M0 = M; mloc0 = mloc; } |
| 568 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
485 | PetscCheck(M==M0,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_INCOMP,"Dimensions of Psplit[%" PetscInt_FMT "] do not match with previous matrices (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,M,M0); |
| 569 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
485 | PetscCheck(mloc==mloc0,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_INCOMP,"Local dimensions of Psplit[%" PetscInt_FMT "] do not match with previous matrices (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,mloc,mloc0); |
| 570 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
485 | PetscCall(PetscObjectReference((PetscObject)Psplit[i])); |
| 571 | } | ||
| 572 | |||
| 573 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
215 | if (st->Psplit) PetscCall(MatDestroyMatrices(st->nsplit,&st->Psplit)); |
| 574 | |||
| 575 | /* allocate space and copy matrices */ | ||
| 576 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
215 | if (n) { |
| 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.
|
215 | PetscCall(PetscMalloc1(n,&st->Psplit)); |
| 578 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
700 | for (i=0;i<n;i++) st->Psplit[i] = Psplit[i]; |
| 579 | } | ||
| 580 | 215 | st->nsplit = n; | |
| 581 | 215 | st->strp = strp; | |
| 582 | 215 | st->state = ST_STATE_INITIAL; | |
| 583 |
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.
|
215 | PetscFunctionReturn(PETSC_SUCCESS); |
| 584 | } | ||
| 585 | |||
| 586 | /*@ | ||
| 587 | STGetSplitPreconditionerTerm - Gets the matrices associated with | ||
| 588 | the split preconditioner. | ||
| 589 | |||
| 590 | Not Collective | ||
| 591 | |||
| 592 | Input Parameters: | ||
| 593 | + st - the spectral transformation context | ||
| 594 | - k - the index of the requested matrix (starting in 0) | ||
| 595 | |||
| 596 | Output Parameter: | ||
| 597 | . Psplit - the returned matrix | ||
| 598 | |||
| 599 | Level: advanced | ||
| 600 | |||
| 601 | .seealso: [](ch:st), `STSetSplitPreconditioner()`, `STGetSplitPreconditionerInfo()` | ||
| 602 | @*/ | ||
| 603 | 145 | PetscErrorCode STGetSplitPreconditionerTerm(ST st,PetscInt k,Mat *Psplit) | |
| 604 | { | ||
| 605 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
145 | PetscFunctionBegin; |
| 606 |
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.
|
145 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 607 |
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.
|
145 | PetscValidLogicalCollectiveInt(st,k,2); |
| 608 |
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.
|
145 | PetscAssertPointer(Psplit,3); |
| 609 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
145 | PetscCheck(k>=0 && k<st->nsplit,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %" PetscInt_FMT,st->nsplit-1); |
| 610 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
145 | PetscCheck(st->Psplit,PetscObjectComm((PetscObject)st),PETSC_ERR_ORDER,"You have not called STSetSplitPreconditioner()"); |
| 611 | 145 | *Psplit = st->Psplit[k]; | |
| 612 |
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.
|
145 | PetscFunctionReturn(PETSC_SUCCESS); |
| 613 | } | ||
| 614 | |||
| 615 | /*@ | ||
| 616 | STGetSplitPreconditionerInfo - Returns the number of matrices of the split | ||
| 617 | preconditioner, as well as the structure flag. | ||
| 618 | |||
| 619 | Not Collective | ||
| 620 | |||
| 621 | Input Parameter: | ||
| 622 | . st - the spectral transformation context | ||
| 623 | |||
| 624 | Output Parameters: | ||
| 625 | + n - the number of matrices passed in `STSetSplitPreconditioner()` | ||
| 626 | - strp - the matrix structure flag passed in `STSetSplitPreconditioner()` | ||
| 627 | |||
| 628 | Level: advanced | ||
| 629 | |||
| 630 | .seealso: [](ch:st), `STSetSplitPreconditioner()`, `STGetSplitPreconditionerTerm()` | ||
| 631 | @*/ | ||
| 632 | 1379 | PetscErrorCode STGetSplitPreconditionerInfo(ST st,PetscInt *n,MatStructure *strp) | |
| 633 | { | ||
| 634 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
1379 | PetscFunctionBegin; |
| 635 |
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.
|
1379 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 636 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 5 times.
|
1379 | if (n) *n = st->nsplit; |
| 637 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1379 | if (strp) *strp = st->strp; |
| 638 |
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.
|
1379 | PetscFunctionReturn(PETSC_SUCCESS); |
| 639 | } | ||
| 640 | |||
| 641 | /*@ | ||
| 642 | STSetShift - Sets the shift associated with the spectral transformation. | ||
| 643 | |||
| 644 | Collective | ||
| 645 | |||
| 646 | Input Parameters: | ||
| 647 | + st - the spectral transformation context | ||
| 648 | - shift - the value of the shift | ||
| 649 | |||
| 650 | Notes: | ||
| 651 | In some spectral transformations, changing the shift may have associated | ||
| 652 | a lot of work, for example recomputing a factorization. | ||
| 653 | |||
| 654 | This function is normally not directly called by users, since the shift is | ||
| 655 | indirectly set by `EPSSetTarget()`. | ||
| 656 | |||
| 657 | Level: intermediate | ||
| 658 | |||
| 659 | .seealso: [](ch:st), `EPSSetTarget()`, `STGetShift()`, `STSetDefaultShift()` | ||
| 660 | @*/ | ||
| 661 | 5990 | PetscErrorCode STSetShift(ST st,PetscScalar shift) | |
| 662 | { | ||
| 663 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
5990 | PetscFunctionBegin; |
| 664 |
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.
|
5990 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 665 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
5990 | PetscValidType(st,1); |
| 666 |
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.
|
5990 | PetscValidLogicalCollectiveScalar(st,shift,2); |
| 667 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
5990 | if (st->sigma != shift) { |
| 668 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
5791 | STCheckNotSeized(st,1); |
| 669 |
7/10✓ 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 taken 8 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
|
5791 | if (st->state==ST_STATE_SETUP) PetscTryTypeMethod(st,setshift,shift); |
| 670 | 5791 | st->sigma = shift; | |
| 671 | } | ||
| 672 | 5990 | st->sigma_set = PETSC_TRUE; | |
| 673 |
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.
|
5990 | PetscFunctionReturn(PETSC_SUCCESS); |
| 674 | } | ||
| 675 | |||
| 676 | /*@ | ||
| 677 | STGetShift - Gets the shift associated with the spectral transformation. | ||
| 678 | |||
| 679 | Not Collective | ||
| 680 | |||
| 681 | Input Parameter: | ||
| 682 | . st - the spectral transformation context | ||
| 683 | |||
| 684 | Output Parameter: | ||
| 685 | . shift - the value of the shift | ||
| 686 | |||
| 687 | Level: intermediate | ||
| 688 | |||
| 689 | .seealso: [](ch:st), `STSetShift()` | ||
| 690 | @*/ | ||
| 691 | 24434 | PetscErrorCode STGetShift(ST st,PetscScalar* shift) | |
| 692 | { | ||
| 693 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
24434 | PetscFunctionBegin; |
| 694 |
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.
|
24434 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 695 |
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.
|
24434 | PetscAssertPointer(shift,2); |
| 696 | 24434 | *shift = st->sigma; | |
| 697 |
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.
|
24434 | PetscFunctionReturn(PETSC_SUCCESS); |
| 698 | } | ||
| 699 | |||
| 700 | /*@ | ||
| 701 | STSetDefaultShift - Sets the value of the shift that should be employed if | ||
| 702 | the user did not specify one. | ||
| 703 | |||
| 704 | Logically Collective | ||
| 705 | |||
| 706 | Input Parameters: | ||
| 707 | + st - the spectral transformation context | ||
| 708 | - defaultshift - the default value of the shift | ||
| 709 | |||
| 710 | Level: developer | ||
| 711 | |||
| 712 | .seealso: [](ch:st), `STSetShift()` | ||
| 713 | @*/ | ||
| 714 | 3303 | PetscErrorCode STSetDefaultShift(ST st,PetscScalar defaultshift) | |
| 715 | { | ||
| 716 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
3303 | PetscFunctionBegin; |
| 717 |
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.
|
3303 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 718 |
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.
|
3303 | PetscValidLogicalCollectiveScalar(st,defaultshift,2); |
| 719 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
3303 | if (st->defsigma != defaultshift) { |
| 720 | 2109 | st->defsigma = defaultshift; | |
| 721 | 2109 | st->state = ST_STATE_INITIAL; | |
| 722 | 2109 | st->opready = PETSC_FALSE; | |
| 723 | } | ||
| 724 |
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.
|
3303 | PetscFunctionReturn(PETSC_SUCCESS); |
| 725 | } | ||
| 726 | |||
| 727 | /*@ | ||
| 728 | STScaleShift - Multiply the shift by a given factor. | ||
| 729 | |||
| 730 | Logically Collective | ||
| 731 | |||
| 732 | Input Parameters: | ||
| 733 | + st - the spectral transformation context | ||
| 734 | - factor - the scaling factor | ||
| 735 | |||
| 736 | Note: | ||
| 737 | This function does not update the transformation matrices, as opposed to | ||
| 738 | `STSetShift()`. | ||
| 739 | |||
| 740 | Level: developer | ||
| 741 | |||
| 742 | .seealso: [](ch:st), `STSetShift()` | ||
| 743 | @*/ | ||
| 744 | 2208 | PetscErrorCode STScaleShift(ST st,PetscScalar factor) | |
| 745 | { | ||
| 746 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2208 | PetscFunctionBegin; |
| 747 |
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.
|
2208 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 748 |
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.
|
2208 | PetscValidLogicalCollectiveScalar(st,factor,2); |
| 749 | 2208 | st->sigma *= factor; | |
| 750 |
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.
|
2208 | PetscFunctionReturn(PETSC_SUCCESS); |
| 751 | } | ||
| 752 | |||
| 753 | /*@ | ||
| 754 | STSetBalanceMatrix - Sets the diagonal matrix to be used for balancing. | ||
| 755 | |||
| 756 | Collective | ||
| 757 | |||
| 758 | Input Parameters: | ||
| 759 | + st - the spectral transformation context | ||
| 760 | - D - the diagonal matrix (represented as a vector) | ||
| 761 | |||
| 762 | Notes: | ||
| 763 | If this matrix is set, `STApply()` will effectively apply $D K^{-1} M D^{-1}$, | ||
| 764 | see discussion at `STGetOperator()`. Use `NULL` to reset a previously passed `D`. | ||
| 765 | |||
| 766 | Balancing is usually set via `EPSSetBalance()`, but the advanced user may use | ||
| 767 | this function to bypass the usual balancing methods. | ||
| 768 | |||
| 769 | Level: developer | ||
| 770 | |||
| 771 | .seealso: [](ch:st), `EPSSetBalance()`, `STApply()`, `STGetBalanceMatrix()`, `STGetOperator()` | ||
| 772 | @*/ | ||
| 773 | 9872 | PetscErrorCode STSetBalanceMatrix(ST st,Vec D) | |
| 774 | { | ||
| 775 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
9872 | PetscFunctionBegin; |
| 776 |
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.
|
9872 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 777 |
8/14✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
|
9872 | if (st->D == D) PetscFunctionReturn(PETSC_SUCCESS); |
| 778 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
130 | STCheckNotSeized(st,1); |
| 779 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
130 | if (D) { |
| 780 |
2/12✗ 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 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
130 | PetscValidHeaderSpecific(D,VEC_CLASSID,2); |
| 781 |
13/32✓ 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 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.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
130 | PetscCheckSameComm(st,1,D,2); |
| 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.
|
130 | PetscCall(PetscObjectReference((PetscObject)D)); |
| 783 | } | ||
| 784 |
4/6✓ Branch 0 taken 2 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(VecDestroy(&st->D)); |
| 785 | 130 | st->D = D; | |
| 786 | 130 | st->state = ST_STATE_INITIAL; | |
| 787 | 130 | st->opready = PETSC_FALSE; | |
| 788 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
130 | PetscFunctionReturn(PETSC_SUCCESS); |
| 789 | } | ||
| 790 | |||
| 791 | /*@ | ||
| 792 | STGetBalanceMatrix - Gets the balance matrix used by the spectral transformation. | ||
| 793 | |||
| 794 | Not Collective | ||
| 795 | |||
| 796 | Input Parameter: | ||
| 797 | . st - the spectral transformation context | ||
| 798 | |||
| 799 | Output Parameter: | ||
| 800 | . D - the diagonal matrix (represented as a vector) | ||
| 801 | |||
| 802 | Note: | ||
| 803 | If the matrix was not set, a `NULL` pointer will be returned. | ||
| 804 | |||
| 805 | Level: developer | ||
| 806 | |||
| 807 | .seealso: [](ch:st), `STSetBalanceMatrix()` | ||
| 808 | @*/ | ||
| 809 | ✗ | PetscErrorCode STGetBalanceMatrix(ST st,Vec *D) | |
| 810 | { | ||
| 811 | ✗ | PetscFunctionBegin; | |
| 812 | ✗ | PetscValidHeaderSpecific(st,ST_CLASSID,1); | |
| 813 | ✗ | PetscAssertPointer(D,2); | |
| 814 | ✗ | *D = st->D; | |
| 815 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
| 816 | } | ||
| 817 | |||
| 818 | /*@ | ||
| 819 | STMatCreateVecs - Get vector(s) compatible with the `ST` matrices. | ||
| 820 | |||
| 821 | Collective | ||
| 822 | |||
| 823 | Input Parameter: | ||
| 824 | . st - the spectral transformation context | ||
| 825 | |||
| 826 | Output Parameters: | ||
| 827 | + right - (optional) vector that the matrix can be multiplied against | ||
| 828 | - left - (optional) vector that the matrix vector product can be stored in | ||
| 829 | |||
| 830 | Level: developer | ||
| 831 | |||
| 832 | .seealso: [](ch:st), `STMatCreateVecsEmpty()` | ||
| 833 | @*/ | ||
| 834 | 4227 | PetscErrorCode STMatCreateVecs(ST st,Vec *right,Vec *left) | |
| 835 | { | ||
| 836 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
4227 | PetscFunctionBegin; |
| 837 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
4227 | STCheckMatrices(st,1); |
| 838 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
4227 | PetscCall(MatCreateVecs(st->A[0],right,left)); |
| 839 |
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.
|
814 | PetscFunctionReturn(PETSC_SUCCESS); |
| 840 | } | ||
| 841 | |||
| 842 | /*@ | ||
| 843 | STMatCreateVecsEmpty - Get vector(s) compatible with the `ST` matrices, i.e., | ||
| 844 | with the same parallel layout, but without internal array. | ||
| 845 | |||
| 846 | Collective | ||
| 847 | |||
| 848 | Input Parameter: | ||
| 849 | . st - the spectral transformation context | ||
| 850 | |||
| 851 | Output Parameters: | ||
| 852 | + right - (optional) vector that the matrix can be multiplied against | ||
| 853 | - left - (optional) vector that the matrix vector product can be stored in | ||
| 854 | |||
| 855 | Level: developer | ||
| 856 | |||
| 857 | .seealso: [](ch:st), `STMatCreateVecs()`, `MatCreateVecsEmpty()` | ||
| 858 | @*/ | ||
| 859 | 9276 | PetscErrorCode STMatCreateVecsEmpty(ST st,Vec *right,Vec *left) | |
| 860 | { | ||
| 861 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
9276 | PetscFunctionBegin; |
| 862 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
9276 | STCheckMatrices(st,1); |
| 863 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9276 | PetscCall(MatCreateVecsEmpty(st->A[0],right,left)); |
| 864 |
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.
|
1812 | PetscFunctionReturn(PETSC_SUCCESS); |
| 865 | } | ||
| 866 | |||
| 867 | /*@ | ||
| 868 | STMatGetSize - Returns the number of rows and columns of the `ST` matrices. | ||
| 869 | |||
| 870 | Not Collective | ||
| 871 | |||
| 872 | Input Parameter: | ||
| 873 | . st - the spectral transformation context | ||
| 874 | |||
| 875 | Output Parameters: | ||
| 876 | + m - the number of global rows | ||
| 877 | - n - the number of global columns | ||
| 878 | |||
| 879 | Level: developer | ||
| 880 | |||
| 881 | .seealso: [](ch:st), `STMatGetLocalSize()` | ||
| 882 | @*/ | ||
| 883 | 9742 | PetscErrorCode STMatGetSize(ST st,PetscInt *m,PetscInt *n) | |
| 884 | { | ||
| 885 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
9742 | PetscFunctionBegin; |
| 886 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
9742 | STCheckMatrices(st,1); |
| 887 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9742 | PetscCall(MatGetSize(st->A[0],m,n)); |
| 888 |
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.
|
1874 | PetscFunctionReturn(PETSC_SUCCESS); |
| 889 | } | ||
| 890 | |||
| 891 | /*@ | ||
| 892 | STMatGetLocalSize - Returns the number of local rows and columns of the `ST` matrices. | ||
| 893 | |||
| 894 | Not Collective | ||
| 895 | |||
| 896 | Input Parameter: | ||
| 897 | . st - the spectral transformation context | ||
| 898 | |||
| 899 | Output Parameters: | ||
| 900 | + m - the number of local rows | ||
| 901 | - n - the number of local columns | ||
| 902 | |||
| 903 | Level: developer | ||
| 904 | |||
| 905 | .seealso: [](ch:st), `STMatGetSize()` | ||
| 906 | @*/ | ||
| 907 | 9742 | PetscErrorCode STMatGetLocalSize(ST st,PetscInt *m,PetscInt *n) | |
| 908 | { | ||
| 909 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
9742 | PetscFunctionBegin; |
| 910 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
9742 | STCheckMatrices(st,1); |
| 911 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
9742 | PetscCall(MatGetLocalSize(st->A[0],m,n)); |
| 912 |
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.
|
1874 | PetscFunctionReturn(PETSC_SUCCESS); |
| 913 | } | ||
| 914 | |||
| 915 | /*@ | ||
| 916 | STSetOptionsPrefix - Sets the prefix used for searching for all | ||
| 917 | `ST` options in the database. | ||
| 918 | |||
| 919 | Logically Collective | ||
| 920 | |||
| 921 | Input Parameters: | ||
| 922 | + st - the spectral transformation context | ||
| 923 | - prefix - the prefix string to prepend to all `ST` option requests | ||
| 924 | |||
| 925 | Notes: | ||
| 926 | A hyphen (-) must NOT be given at the beginning of the prefix name. | ||
| 927 | The first character of all runtime options is AUTOMATICALLY the | ||
| 928 | hyphen. | ||
| 929 | |||
| 930 | Level: advanced | ||
| 931 | |||
| 932 | .seealso: [](ch:st), `STAppendOptionsPrefix()`, `STGetOptionsPrefix()` | ||
| 933 | @*/ | ||
| 934 | 2413 | PetscErrorCode STSetOptionsPrefix(ST st,const char prefix[]) | |
| 935 | { | ||
| 936 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2413 | PetscFunctionBegin; |
| 937 |
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.
|
2413 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 938 |
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.
|
2413 | if (!st->ksp) PetscCall(STGetKSP(st,&st->ksp)); |
| 939 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2413 | PetscCall(KSPSetOptionsPrefix(st->ksp,prefix)); |
| 940 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2413 | PetscCall(KSPAppendOptionsPrefix(st->ksp,"st_")); |
| 941 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2413 | PetscCall(PetscObjectSetOptionsPrefix((PetscObject)st,prefix)); |
| 942 |
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.
|
464 | PetscFunctionReturn(PETSC_SUCCESS); |
| 943 | } | ||
| 944 | |||
| 945 | /*@ | ||
| 946 | STAppendOptionsPrefix - Appends to the prefix used for searching for all | ||
| 947 | `ST` options in the database. | ||
| 948 | |||
| 949 | Logically Collective | ||
| 950 | |||
| 951 | Input Parameters: | ||
| 952 | + st - the spectral transformation context | ||
| 953 | - prefix - the prefix string to prepend to all `ST` option requests | ||
| 954 | |||
| 955 | Notes: | ||
| 956 | A hyphen (-) must NOT be given at the beginning of the prefix name. | ||
| 957 | The first character of all runtime options is AUTOMATICALLY the | ||
| 958 | hyphen. | ||
| 959 | |||
| 960 | Level: advanced | ||
| 961 | |||
| 962 | .seealso: [](ch:st), `STSetOptionsPrefix()`, `STGetOptionsPrefix()` | ||
| 963 | @*/ | ||
| 964 | 1934 | PetscErrorCode STAppendOptionsPrefix(ST st,const char prefix[]) | |
| 965 | { | ||
| 966 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
1934 | PetscFunctionBegin; |
| 967 |
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.
|
1934 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 968 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1934 | PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)st,prefix)); |
| 969 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
1934 | if (!st->ksp) PetscCall(STGetKSP(st,&st->ksp)); |
| 970 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1934 | PetscCall(KSPSetOptionsPrefix(st->ksp,((PetscObject)st)->prefix)); |
| 971 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
1934 | PetscCall(KSPAppendOptionsPrefix(st->ksp,"st_")); |
| 972 |
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.
|
382 | PetscFunctionReturn(PETSC_SUCCESS); |
| 973 | } | ||
| 974 | |||
| 975 | /*@ | ||
| 976 | STGetOptionsPrefix - Gets the prefix used for searching for all | ||
| 977 | `ST` options in the database. | ||
| 978 | |||
| 979 | Not Collective | ||
| 980 | |||
| 981 | Input Parameter: | ||
| 982 | . st - the spectral transformation context | ||
| 983 | |||
| 984 | Output Parameter: | ||
| 985 | . prefix - pointer to the prefix string used, is returned | ||
| 986 | |||
| 987 | Level: advanced | ||
| 988 | |||
| 989 | .seealso: [](ch:st), `STSetOptionsPrefix()`, `STAppendOptionsPrefix()` | ||
| 990 | @*/ | ||
| 991 | 21 | PetscErrorCode STGetOptionsPrefix(ST st,const char *prefix[]) | |
| 992 | { | ||
| 993 |
0/2✗ Branch 0 not taken.
✗ Branch 1 not taken.
|
21 | PetscFunctionBegin; |
| 994 |
0/16✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ 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.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
21 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 995 |
0/8✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
21 | PetscAssertPointer(prefix,2); |
| 996 |
1/6✗ 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.
|
21 | PetscCall(PetscObjectGetOptionsPrefix((PetscObject)st,prefix)); |
| 997 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
| 998 | } | ||
| 999 | |||
| 1000 | /*@ | ||
| 1001 | STView - Prints the `ST` data structure. | ||
| 1002 | |||
| 1003 | Collective | ||
| 1004 | |||
| 1005 | Input Parameters: | ||
| 1006 | + st - the ST context | ||
| 1007 | - viewer - optional visualization context | ||
| 1008 | |||
| 1009 | Note: | ||
| 1010 | The available visualization contexts include | ||
| 1011 | + `PETSC_VIEWER_STDOUT_SELF` - standard output (default) | ||
| 1012 | - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the | ||
| 1013 | first process opens the file; all other processes send their data to the | ||
| 1014 | first one to print | ||
| 1015 | |||
| 1016 | The user can open an alternative visualization context with `PetscViewerASCIIOpen()` | ||
| 1017 | to output to a specified file. | ||
| 1018 | |||
| 1019 | Use `STViewFromOptions()` to allow the user to select many different `PetscViewerType` | ||
| 1020 | and formats from the options database. | ||
| 1021 | |||
| 1022 | Level: beginner | ||
| 1023 | |||
| 1024 | .seealso: [](ch:st), `STCreate()`, `STViewFromOptions()` | ||
| 1025 | @*/ | ||
| 1026 | 126 | PetscErrorCode STView(ST st,PetscViewer viewer) | |
| 1027 | { | ||
| 1028 | 126 | STType cstr; | |
| 1029 | 126 | char str[50]; | |
| 1030 | 126 | PetscBool isascii,isstring; | |
| 1031 | |||
| 1032 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
126 | PetscFunctionBegin; |
| 1033 |
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.
|
126 | PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| 1034 |
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.
|
126 | if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)st),&viewer)); |
| 1035 |
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.
|
126 | PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); |
| 1036 |
13/32✓ 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 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.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
126 | PetscCheckSameComm(st,1,viewer,2); |
| 1037 | |||
| 1038 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii)); |
| 1039 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring)); |
| 1040 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
126 | if (isascii) { |
| 1041 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)st,viewer)); |
| 1042 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(PetscViewerASCIIPushTab(viewer)); |
| 1043 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
126 | PetscTryTypeMethod(st,view,viewer); |
| 1044 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(PetscViewerASCIIPopTab(viewer)); |
| 1045 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(SlepcSNPrintfScalar(str,sizeof(str),st->sigma,PETSC_FALSE)); |
| 1046 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(PetscViewerASCIIPrintf(viewer," shift: %s\n",str)); |
| 1047 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
126 | PetscCall(PetscViewerASCIIPrintf(viewer," number of matrices: %" PetscInt_FMT "\n",st->nmat)); |
| 1048 |
2/3✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
|
126 | switch (st->matmode) { |
| 1049 | case ST_MATMODE_COPY: | ||
| 1050 | break; | ||
| 1051 | ✗ | case ST_MATMODE_INPLACE: | |
| 1052 | ✗ | PetscCall(PetscViewerASCIIPrintf(viewer," shifting the matrix and unshifting at exit\n")); | |
| 1053 | break; | ||
| 1054 | 10 | case ST_MATMODE_SHELL: | |
| 1055 |
4/6✓ Branch 0 taken 2 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(PetscViewerASCIIPrintf(viewer," using a shell matrix\n")); |
| 1056 | break; | ||
| 1057 | } | ||
| 1058 |
7/10✓ 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 taken 8 times.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
|
126 | if (st->nmat>1 && st->matmode != ST_MATMODE_SHELL) PetscCall(PetscViewerASCIIPrintf(viewer," nonzero pattern of the matrices: %s\n",MatStructures[st->str])); |
| 1059 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
126 | if (st->Psplit) PetscCall(PetscViewerASCIIPrintf(viewer," using split preconditioner matrices with %s\n",MatStructures[st->strp])); |
| 1060 |
3/10✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
|
126 | if (st->transform && st->nmat>2) PetscCall(PetscViewerASCIIPrintf(viewer," computing transformed matrices\n")); |
| 1061 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
126 | if (st->structured) PetscCall(PetscViewerASCIIPrintf(viewer," exploiting structure in the application of the operator\n")); |
| 1062 | ✗ | } else if (isstring) { | |
| 1063 | ✗ | PetscCall(STGetType(st,&cstr)); | |
| 1064 | ✗ | PetscCall(PetscViewerStringSPrintf(viewer," %-7.7s",cstr)); | |
| 1065 | ✗ | PetscTryTypeMethod(st,view,viewer); | |
| 1066 | } | ||
| 1067 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
126 | if (st->usesksp) { |
| 1068 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
59 | if (!st->ksp) PetscCall(STGetKSP(st,&st->ksp)); |
| 1069 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
59 | PetscCall(PetscViewerASCIIPushTab(viewer)); |
| 1070 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
59 | PetscCall(KSPView(st->ksp,viewer)); |
| 1071 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
59 | PetscCall(PetscViewerASCIIPopTab(viewer)); |
| 1072 | } | ||
| 1073 |
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.
|
24 | PetscFunctionReturn(PETSC_SUCCESS); |
| 1074 | } | ||
| 1075 | |||
| 1076 | /*@ | ||
| 1077 | STViewFromOptions - View (print) an `ST` object based on values in the options database | ||
| 1078 | |||
| 1079 | Collective | ||
| 1080 | |||
| 1081 | Input Parameters: | ||
| 1082 | + st - the spectral transformation context | ||
| 1083 | . obj - optional object that provides the options prefix used to query the options database | ||
| 1084 | - name - command line option | ||
| 1085 | |||
| 1086 | Level: intermediate | ||
| 1087 | |||
| 1088 | .seealso: [](ch:st), `STView()`, `STCreate()`, `PetscObjectViewFromOptions()` | ||
| 1089 | @*/ | ||
| 1090 | ✗ | PetscErrorCode STViewFromOptions(ST st,PetscObject obj,const char name[]) | |
| 1091 | { | ||
| 1092 | ✗ | PetscFunctionBegin; | |
| 1093 | ✗ | PetscValidHeaderSpecific(st,ST_CLASSID,1); | |
| 1094 | ✗ | PetscCall(PetscObjectViewFromOptions((PetscObject)st,obj,name)); | |
| 1095 | ✗ | PetscFunctionReturn(PETSC_SUCCESS); | |
| 1096 | } | ||
| 1097 | |||
| 1098 | /*@C | ||
| 1099 | STRegister - Adds a method to the spectral transformation package. | ||
| 1100 | |||
| 1101 | Not Collective | ||
| 1102 | |||
| 1103 | Input Parameters: | ||
| 1104 | + name - name of a new user-defined transformation | ||
| 1105 | - function - routine to create method | ||
| 1106 | |||
| 1107 | Notes: | ||
| 1108 | `STRegister()` may be called multiple times to add several user-defined | ||
| 1109 | spectral transformations. | ||
| 1110 | |||
| 1111 | Example Usage: | ||
| 1112 | .vb | ||
| 1113 | STRegister("my_transform",MyTransformCreate); | ||
| 1114 | .ve | ||
| 1115 | |||
| 1116 | Then, your spectral transform can be chosen with the procedural interface via | ||
| 1117 | .vb | ||
| 1118 | STSetType(st,"my_transform") | ||
| 1119 | .ve | ||
| 1120 | or at runtime via the option `-st_type my_transform` | ||
| 1121 | |||
| 1122 | Level: advanced | ||
| 1123 | |||
| 1124 | .seealso: [](ch:st), `STSetType()`, `STRegisterAll()` | ||
| 1125 | @*/ | ||
| 1126 | 47994 | PetscErrorCode STRegister(const char *name,PetscErrorCode (*function)(ST)) | |
| 1127 | { | ||
| 1128 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
47994 | PetscFunctionBegin; |
| 1129 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
47994 | PetscCall(STInitializePackage()); |
| 1130 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
47994 | PetscCall(PetscFunctionListAdd(&STList,name,function)); |
| 1131 |
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.
|
9426 | PetscFunctionReturn(PETSC_SUCCESS); |
| 1132 | } | ||
| 1133 |