| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | /* | ||
| 2 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 3 | SLEPc - Scalable Library for Eigenvalue Problem Computations | ||
| 4 | Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain | ||
| 5 | |||
| 6 | This file is part of SLEPc. | ||
| 7 | SLEPc is distributed under a 2-clause BSD license (see LICENSE). | ||
| 8 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 9 | */ | ||
| 10 | |||
| 11 | static char help[] = "Linear Response eigenvalue problem.\n\n" | ||
| 12 | "The command line options are:\n" | ||
| 13 | " -n <n>, where <n> = dimension of the blocks.\n" | ||
| 14 | " -reduced <0/1>, to use the reduced form of the LREP.\n" | ||
| 15 | " -A <filename>, A matrix.\n" | ||
| 16 | " -B <filename>, B matrix.\n\n"; | ||
| 17 | |||
| 18 | #include <slepceps.h> | ||
| 19 | |||
| 20 | /* | ||
| 21 | This example is similar to ex55.c, but with a real matrix | ||
| 22 | |||
| 23 | H = [ A B | ||
| 24 | -B -A ], | ||
| 25 | |||
| 26 | where A,B are real symmetric with the following Toeplitz structure: | ||
| 27 | |||
| 28 | A = pentadiag{a,b,c,b,a} | ||
| 29 | B = tridiag{b,d,b} | ||
| 30 | |||
| 31 | where a,b,c,d are real values. | ||
| 32 | |||
| 33 | Alternatively, use -A and -B command-line options to load matrices from file. | ||
| 34 | */ | ||
| 35 | |||
| 36 | 154 | int main(int argc,char **argv) | |
| 37 | { | ||
| 38 | 154 | Mat H,A,B,K,M; /* problem matrices */ | |
| 39 | 154 | EPS eps; /* eigenproblem solver context */ | |
| 40 | 154 | PetscScalar a,b,c,d; | |
| 41 | 154 | PetscReal lev; | |
| 42 | 154 | PetscInt n=24,Istart,Iend,i,nconv; | |
| 43 | 154 | PetscBool flg,terse,checkorthog,nest=PETSC_FALSE,reduced=PETSC_FALSE; | |
| 44 | 154 | Vec t,*x,*y; | |
| 45 | 154 | PetscViewer viewer; | |
| 46 | 154 | char filename[PETSC_MAX_PATH_LEN]; | |
| 47 | |||
| 48 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
154 | PetscFunctionBeginUser; |
| 49 |
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.
|
154 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
| 50 | 154 | PetscCheck(!PetscDefined(USE_COMPLEX),PETSC_COMM_SELF,PETSC_ERR_SUP,"This example cannot be used with complex scalars"); | |
| 51 | |||
| 52 |
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.
|
154 | PetscCall(PetscOptionsGetBool(NULL,NULL,"-reduced",&reduced,NULL)); |
| 53 | |||
| 54 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 55 | Compute (or load) the problem matrices A and B | ||
| 56 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 57 | |||
| 58 |
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.
|
154 | PetscCall(PetscOptionsHasName(NULL,NULL,"-A",&flg)); |
| 59 | |||
| 60 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
|
154 | if (flg) { |
| 61 | |||
| 62 | ✗ | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nLinear Response Eigenvalue Problem stored in file%s\n\n",reduced?" (reduced)":"")); | |
| 63 | ✗ | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrices from a binary file...\n")); | |
| 64 | |||
| 65 | ✗ | PetscCall(PetscOptionsGetString(NULL,NULL,"-A",filename,sizeof(filename),&flg)); | |
| 66 | ✗ | PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer)); | |
| 67 | ✗ | PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); | |
| 68 | ✗ | PetscCall(MatSetFromOptions(A)); | |
| 69 | ✗ | PetscCall(MatLoad(A,viewer)); | |
| 70 | ✗ | PetscCall(PetscViewerDestroy(&viewer)); | |
| 71 | |||
| 72 | ✗ | PetscCall(PetscOptionsGetString(NULL,NULL,"-B",filename,sizeof(filename),&flg)); | |
| 73 | ✗ | PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Must indicate file for both A and B matrices"); | |
| 74 | ✗ | PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer)); | |
| 75 | ✗ | PetscCall(MatCreate(PETSC_COMM_WORLD,&B)); | |
| 76 | ✗ | PetscCall(MatSetFromOptions(B)); | |
| 77 | ✗ | PetscCall(MatLoad(B,viewer)); | |
| 78 | ✗ | PetscCall(PetscViewerDestroy(&viewer)); | |
| 79 | |||
| 80 | } else { | ||
| 81 | |||
| 82 |
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.
|
154 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); |
| 83 |
7/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✓ Branch 5 taken 1 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
231 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nLinear Response Eigenvalue Problem, n=%" PetscInt_FMT "%s\n\n",n,reduced?" (reduced)":"")); |
| 84 | |||
| 85 | 154 | a = -0.1; | |
| 86 | 154 | b = 1.0; | |
| 87 | 154 | c = 4.5; | |
| 88 | 154 | d = 2.0; | |
| 89 | |||
| 90 |
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.
|
154 | PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); |
| 91 |
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.
|
154 | PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n)); |
| 92 |
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.
|
154 | PetscCall(MatSetFromOptions(A)); |
| 93 | |||
| 94 |
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.
|
154 | PetscCall(MatCreate(PETSC_COMM_WORLD,&B)); |
| 95 |
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.
|
154 | PetscCall(MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,n,n)); |
| 96 |
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.
|
154 | PetscCall(MatSetFromOptions(B)); |
| 97 | |||
| 98 |
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.
|
154 | PetscCall(MatGetOwnershipRange(A,&Istart,&Iend)); |
| 99 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
3358 | for (i=Istart;i<Iend;i++) { |
| 100 |
6/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
3204 | if (i>1) PetscCall(MatSetValue(A,i,i-2,a,INSERT_VALUES)); |
| 101 |
6/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
3204 | if (i>0) PetscCall(MatSetValue(A,i,i-1,b,INSERT_VALUES)); |
| 102 |
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.
|
3204 | PetscCall(MatSetValue(A,i,i,c,INSERT_VALUES)); |
| 103 |
6/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
3204 | if (i<n-1) PetscCall(MatSetValue(A,i,i+1,b,INSERT_VALUES)); |
| 104 |
6/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
3204 | if (i<n-2) PetscCall(MatSetValue(A,i,i+2,a,INSERT_VALUES)); |
| 105 | } | ||
| 106 | |||
| 107 |
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.
|
154 | PetscCall(MatGetOwnershipRange(B,&Istart,&Iend)); |
| 108 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
3358 | for (i=Istart;i<Iend;i++) { |
| 109 |
6/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
3204 | if (i>0) PetscCall(MatSetValue(B,i,i-1,b,INSERT_VALUES)); |
| 110 |
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.
|
3204 | PetscCall(MatSetValue(B,i,i,d,INSERT_VALUES)); |
| 111 |
6/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
3204 | if (i<n-1) PetscCall(MatSetValue(B,i,i+1,b,INSERT_VALUES)); |
| 112 | } | ||
| 113 | |||
| 114 |
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.
|
154 | PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); |
| 115 |
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.
|
154 | PetscCall(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY)); |
| 116 |
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.
|
154 | PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); |
| 117 |
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.
|
154 | PetscCall(MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY)); |
| 118 | } | ||
| 119 | |||
| 120 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
154 | if (reduced) { |
| 121 | |||
| 122 | 77 | PetscInt ma,na,Ma,Na; | |
| 123 | 77 | const PetscScalar scal[] = { 1.0, -1.0 }; | |
| 124 | |||
| 125 |
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.
|
77 | PetscCall(MatGetSize(A,&Ma,&Na)); |
| 126 |
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.
|
77 | PetscCall(MatGetLocalSize(A,&ma,&na)); |
| 127 | |||
| 128 | /* K = A-B */ | ||
| 129 |
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.
|
77 | PetscCall(MatCreate(PetscObjectComm((PetscObject)A),&K)); |
| 130 |
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.
|
77 | PetscCall(MatSetSizes(K,ma,na,Ma,Na)); |
| 131 |
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.
|
77 | PetscCall(MatSetType(K,MATCOMPOSITE)); |
| 132 |
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.
|
77 | PetscCall(MatCompositeAddMat(K,A)); |
| 133 |
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.
|
77 | PetscCall(MatCompositeAddMat(K,B)); |
| 134 |
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.
|
77 | PetscCall(MatAssemblyBegin(K,MAT_FINAL_ASSEMBLY)); |
| 135 |
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.
|
77 | PetscCall(MatAssemblyEnd(K,MAT_FINAL_ASSEMBLY)); |
| 136 |
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.
|
77 | PetscCall(MatCompositeSetScalings(K,scal)); |
| 137 | |||
| 138 | /* M = A+B */ | ||
| 139 |
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.
|
77 | PetscCall(MatCreate(PetscObjectComm((PetscObject)A),&M)); |
| 140 |
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.
|
77 | PetscCall(MatSetSizes(M,ma,na,Ma,Na)); |
| 141 |
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.
|
77 | PetscCall(MatSetType(M,MATCOMPOSITE)); |
| 142 |
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.
|
77 | PetscCall(MatCompositeAddMat(M,A)); |
| 143 |
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.
|
77 | PetscCall(MatCompositeAddMat(M,B)); |
| 144 |
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.
|
77 | PetscCall(MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY)); |
| 145 |
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.
|
77 | PetscCall(MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY)); |
| 146 | |||
| 147 |
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.
|
77 | PetscCall(MatCreateLREP(K,M,reduced,&H)); |
| 148 |
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.
|
77 | PetscCall(MatDestroy(&K)); |
| 149 |
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.
|
77 | PetscCall(MatDestroy(&M)); |
| 150 | |||
| 151 |
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.
|
77 | } else PetscCall(MatCreateLREP(A,B,reduced,&H)); |
| 152 | |||
| 153 | /* if you prefer, set the vector type so that MatCreateVecs() returns nested vectors */ | ||
| 154 |
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.
|
154 | PetscCall(PetscOptionsGetBool(NULL,NULL,"-nest",&nest,NULL)); |
| 155 |
6/8✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
154 | if (nest) PetscCall(MatNestSetVecType(H,VECNEST)); |
| 156 | |||
| 157 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 158 | Create the eigensolver and set various options | ||
| 159 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 160 | |||
| 161 |
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.
|
154 | PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps)); |
| 162 |
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.
|
154 | PetscCall(EPSSetOperators(eps,H,NULL)); |
| 163 |
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.
|
154 | PetscCall(EPSSetProblemType(eps,EPS_LREP)); |
| 164 |
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.
|
154 | PetscCall(EPSSetFromOptions(eps)); |
| 165 | |||
| 166 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 167 | Solve the eigensystem and display solution | ||
| 168 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 169 | |||
| 170 |
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.
|
154 | PetscCall(EPSSolve(eps)); |
| 171 | |||
| 172 | /* show detailed info unless -terse option is given by user */ | ||
| 173 |
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.
|
154 | PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
| 174 |
5/8✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
154 | if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL)); |
| 175 | else { | ||
| 176 | ✗ | PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL)); | |
| 177 | ✗ | PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD)); | |
| 178 | ✗ | PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD)); | |
| 179 | ✗ | PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); | |
| 180 | } | ||
| 181 | |||
| 182 | /* check bi-orthogonality */ | ||
| 183 |
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.
|
154 | PetscCall(PetscOptionsHasName(NULL,NULL,"-checkorthog",&checkorthog)); |
| 184 |
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.
|
154 | PetscCall(EPSGetConverged(eps,&nconv)); |
| 185 |
2/4✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 5 times.
✗ Branch 3 not taken.
|
154 | if (checkorthog && nconv>0) { |
| 186 |
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.
|
154 | PetscCall(MatCreateVecs(H,&t,NULL)); |
| 187 |
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.
|
154 | PetscCall(VecDuplicateVecs(t,nconv,&x)); |
| 188 |
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.
|
154 | PetscCall(VecDuplicateVecs(t,nconv,&y)); |
| 189 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
2330 | for (i=0;i<nconv;i++) { |
| 190 |
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.
|
2176 | PetscCall(EPSGetEigenvector(eps,i,x[i],NULL)); |
| 191 |
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.
|
2176 | PetscCall(EPSGetLeftEigenvector(eps,i,y[i],NULL)); |
| 192 | } | ||
| 193 |
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.
|
154 | PetscCall(VecCheckOrthogonality(x,nconv,y,nconv,NULL,NULL,&lev)); |
| 194 |
5/8✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
154 | if (lev<100*PETSC_MACHINE_EPSILON) PetscCall(PetscPrintf(PETSC_COMM_WORLD," Level of bi-orthogonality of eigenvectors < 100*eps\n\n")); |
| 195 | ✗ | else PetscCall(PetscPrintf(PETSC_COMM_WORLD," Level of bi-orthogonality of eigenvectors: %g\n\n",(double)lev)); | |
| 196 |
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.
|
154 | PetscCall(VecDestroy(&t)); |
| 197 |
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.
|
154 | PetscCall(VecDestroyVecs(nconv,&x)); |
| 198 |
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.
|
154 | PetscCall(VecDestroyVecs(nconv,&y)); |
| 199 | } | ||
| 200 | |||
| 201 |
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.
|
154 | PetscCall(EPSDestroy(&eps)); |
| 202 |
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.
|
154 | PetscCall(MatDestroy(&A)); |
| 203 |
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.
|
154 | PetscCall(MatDestroy(&B)); |
| 204 |
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.
|
154 | PetscCall(MatDestroy(&H)); |
| 205 |
2/6✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
154 | PetscCall(SlepcFinalize()); |
| 206 | return 0; | ||
| 207 | } | ||
| 208 | |||
| 209 | /*TEST | ||
| 210 | |||
| 211 | build: | ||
| 212 | requires: !complex | ||
| 213 | |||
| 214 | testset: | ||
| 215 | args: -eps_nev 4 -eps_ncv 16 -terse -checkorthog -reduced {{0 1}} -nest {{0 1}} | ||
| 216 | nsize: {{1 2}} | ||
| 217 | filter: sed -e "s/ (reduced)//" | ||
| 218 | output_file: output/ex58_1.out | ||
| 219 | test: | ||
| 220 | suffix: 1 | ||
| 221 | test: | ||
| 222 | suffix: 1_dense | ||
| 223 | args: -mat_type dense | ||
| 224 | test: | ||
| 225 | suffix: 1_cuda | ||
| 226 | args: -mat_type aijcusparse | ||
| 227 | requires: cuda | ||
| 228 | test: | ||
| 229 | suffix: 1_hip | ||
| 230 | args: -mat_type aijhipsparse | ||
| 231 | requires: hip | ||
| 232 | |||
| 233 | test: | ||
| 234 | args: -n 90 -eps_threshold_absolute 2.4 -eps_ncv 10 -terse -checkorthog -reduced {{0 1}} | ||
| 235 | filter: sed -e "s/ (reduced)//" | ||
| 236 | suffix: 2 | ||
| 237 | |||
| 238 | TEST*/ | ||
| 239 |