| 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[] = "User-defined split preconditioner when solving a generalized eigenproblem.\n\n" | ||
| 12 | "The command line options are:\n" | ||
| 13 | " -n <n>, where <n> = number of grid subdivisions in x dimension.\n" | ||
| 14 | " -m <m>, where <m> = number of grid subdivisions in y dimension.\n\n"; | ||
| 15 | |||
| 16 | #include <slepceps.h> | ||
| 17 | |||
| 18 | 70 | int main(int argc,char **argv) | |
| 19 | { | ||
| 20 | 70 | Mat A,B,A0,B0,mats[2]; /* problem matrices and sparser approximations */ | |
| 21 | 70 | EPS eps; /* eigenproblem solver context */ | |
| 22 | 70 | ST st; | |
| 23 | 70 | PetscInt N,n=24,m,Istart,Iend,II,i,j; | |
| 24 | 70 | PetscBool flag,terse; | |
| 25 | |||
| 26 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
70 | PetscFunctionBeginUser; |
| 27 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
| 28 | |||
| 29 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); |
| 30 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag)); |
| 31 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
70 | if (!flag) m=n; |
| 32 | 70 | N = n*m; | |
| 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.
|
70 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nGHEP with split preconditioner, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid)\n\n",N,n,m)); |
| 34 | |||
| 35 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 36 | Compute the problem matrices A and B | ||
| 37 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 38 | |||
| 39 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); |
| 40 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
| 41 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatSetFromOptions(A)); |
| 42 | |||
| 43 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatCreate(PETSC_COMM_WORLD,&B)); |
| 44 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
| 45 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatSetFromOptions(B)); |
| 46 | |||
| 47 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatGetOwnershipRange(A,&Istart,&Iend)); |
| 48 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
34630 | for (II=Istart;II<Iend;II++) { |
| 49 | 34560 | i = II/n; j = II-i*n; | |
| 50 |
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.
|
34560 | if (i>0) PetscCall(MatSetValue(A,II,II-n,-0.2,INSERT_VALUES)); |
| 51 |
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.
|
34560 | if (i<m-1) PetscCall(MatSetValue(A,II,II+n,-0.2,INSERT_VALUES)); |
| 52 |
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.
|
34560 | if (j>0) PetscCall(MatSetValue(A,II,II-1,-3.0,INSERT_VALUES)); |
| 53 |
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.
|
34560 | if (j<n-1) PetscCall(MatSetValue(A,II,II+1,-3.0,INSERT_VALUES)); |
| 54 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
34560 | PetscCall(MatSetValue(A,II,II,7.0,INSERT_VALUES)); |
| 55 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
34560 | PetscCall(MatSetValue(B,II,II,2.0,INSERT_VALUES)); |
| 56 | } | ||
| 57 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
70 | if (Istart==0) { |
| 58 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatSetValue(B,0,0,6.0,INSERT_VALUES)); |
| 59 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatSetValue(B,0,1,-1.0,INSERT_VALUES)); |
| 60 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatSetValue(B,1,0,-1.0,INSERT_VALUES)); |
| 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.
|
60 | PetscCall(MatSetValue(B,1,1,1.0,INSERT_VALUES)); |
| 62 | } | ||
| 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.
|
70 | PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); |
| 64 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); |
| 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.
|
70 | PetscCall(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY)); |
| 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.
|
70 | PetscCall(MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY)); |
| 67 | |||
| 68 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 69 | Compute sparser approximations A0 and B0 | ||
| 70 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 71 | |||
| 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.
|
70 | PetscCall(MatCreate(PETSC_COMM_WORLD,&A0)); |
| 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.
|
70 | PetscCall(MatSetSizes(A0,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
| 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.
|
70 | PetscCall(MatSetFromOptions(A0)); |
| 75 | |||
| 76 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatCreate(PETSC_COMM_WORLD,&B0)); |
| 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.
|
70 | PetscCall(MatSetSizes(B0,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
| 78 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatSetFromOptions(B0)); |
| 79 | |||
| 80 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatGetOwnershipRange(A0,&Istart,&Iend)); |
| 81 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
34630 | for (II=Istart;II<Iend;II++) { |
| 82 | 34560 | i = II/n; j = II-i*n; | |
| 83 |
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.
|
34560 | if (j>0) PetscCall(MatSetValue(A0,II,II-1,-3.0,INSERT_VALUES)); |
| 84 |
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.
|
34560 | if (j<n-1) PetscCall(MatSetValue(A0,II,II+1,-3.0,INSERT_VALUES)); |
| 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.
|
34560 | PetscCall(MatSetValue(A0,II,II,7.0,INSERT_VALUES)); |
| 86 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
34560 | PetscCall(MatSetValue(B0,II,II,2.0,INSERT_VALUES)); |
| 87 | } | ||
| 88 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
70 | if (Istart==0) { |
| 89 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatSetValue(B0,0,0,6.0,INSERT_VALUES)); |
| 90 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatSetValue(B0,1,1,1.0,INSERT_VALUES)); |
| 91 | } | ||
| 92 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatAssemblyBegin(A0,MAT_FINAL_ASSEMBLY)); |
| 93 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatAssemblyEnd(A0,MAT_FINAL_ASSEMBLY)); |
| 94 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatAssemblyBegin(B0,MAT_FINAL_ASSEMBLY)); |
| 95 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatAssemblyEnd(B0,MAT_FINAL_ASSEMBLY)); |
| 96 | |||
| 97 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 98 | Create the eigensolver and set various options | ||
| 99 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 100 | |||
| 101 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps)); |
| 102 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(EPSSetOperators(eps,A,B)); |
| 103 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(EPSSetProblemType(eps,EPS_GHEP)); |
| 104 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(EPSGetST(eps,&st)); |
| 105 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(STSetType(st,STSINVERT)); |
| 106 | 70 | mats[0] = A0; mats[1] = B0; | |
| 107 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(STSetSplitPreconditioner(st,2,mats,SUBSET_NONZERO_PATTERN)); |
| 108 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(EPSSetTarget(eps,0.0)); |
| 109 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE)); |
| 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.
|
70 | PetscCall(EPSSetFromOptions(eps)); |
| 111 | |||
| 112 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 113 | Solve the eigensystem and display solution | ||
| 114 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 115 | |||
| 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.
|
70 | PetscCall(EPSSolve(eps)); |
| 117 | |||
| 118 | /* show detailed info unless -terse option is given by user */ | ||
| 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.
|
70 | PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
| 120 |
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.
|
70 | if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL)); |
| 121 | else { | ||
| 122 | ✗ | PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL)); | |
| 123 | ✗ | PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD)); | |
| 124 | ✗ | PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD)); | |
| 125 | ✗ | PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); | |
| 126 | } | ||
| 127 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(EPSDestroy(&eps)); |
| 128 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatDestroy(&A)); |
| 129 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatDestroy(&B)); |
| 130 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatDestroy(&A0)); |
| 131 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
70 | PetscCall(MatDestroy(&B0)); |
| 132 |
3/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
70 | PetscCall(SlepcFinalize()); |
| 133 | return 0; | ||
| 134 | } | ||
| 135 | |||
| 136 | /*TEST | ||
| 137 | |||
| 138 | testset: | ||
| 139 | args: -eps_nev 4 -terse | ||
| 140 | output_file: output/ex49_1.out | ||
| 141 | requires: !single | ||
| 142 | test: | ||
| 143 | suffix: 1 | ||
| 144 | test: | ||
| 145 | suffix: 1_jd | ||
| 146 | args: -eps_type jd -st_type precond | ||
| 147 | test: | ||
| 148 | suffix: 1_lobpcg | ||
| 149 | args: -eps_type lobpcg -st_type precond -eps_smallest_real -st_shift 0.2 | ||
| 150 | |||
| 151 | testset: | ||
| 152 | args: -eps_type ciss -eps_all -rg_type ellipse -rg_ellipse_center 0 -rg_ellipse_radius 0.34 -rg_ellipse_vscale .2 -st_ksp_type gmres -terse | ||
| 153 | output_file: output/ex49_2.out | ||
| 154 | test: | ||
| 155 | suffix: 2 | ||
| 156 | test: | ||
| 157 | suffix: 2_nost | ||
| 158 | args: -eps_ciss_usest 0 | ||
| 159 | requires: !single | ||
| 160 | test: | ||
| 161 | suffix: 2_par | ||
| 162 | nsize: 2 | ||
| 163 | args: -eps_ciss_partitions 2 | ||
| 164 | requires: !single | ||
| 165 | |||
| 166 | TEST*/ | ||
| 167 |