| 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[] = "Illustrates the use of a region for filtering; the number of wanted eigenvalues is not known a priori.\n\n" | ||
| 12 | "The problem is the Brusselator wave model as in ex9.c.\n" | ||
| 13 | "The command line options are:\n" | ||
| 14 | " -n <n>, where <n> = block dimension of the 2x2 block matrix.\n" | ||
| 15 | " -L <L>, where <L> = bifurcation parameter.\n" | ||
| 16 | " -alpha <alpha>, -beta <beta>, -delta1 <delta1>, -delta2 <delta2>,\n" | ||
| 17 | " where <alpha> <beta> <delta1> <delta2> = model parameters.\n\n"; | ||
| 18 | |||
| 19 | #include <slepceps.h> | ||
| 20 | |||
| 21 | /* | ||
| 22 | This example tries to compute all eigenvalues lying outside the real axis. | ||
| 23 | This could be achieved by computing LARGEST_IMAGINARY eigenvalues, but | ||
| 24 | here we take a different route: define a region of the complex plane where | ||
| 25 | eigenvalues must be emphasized (eigenvalues outside the region are filtered | ||
| 26 | out). In this case, we select the region as the complement of a thin stripe | ||
| 27 | around the real axis. | ||
| 28 | */ | ||
| 29 | |||
| 30 | PetscErrorCode MatMult_Brussel(Mat,Vec,Vec); | ||
| 31 | PetscErrorCode MyStoppingTest(EPS,PetscInt,PetscInt,PetscInt,PetscInt,EPSConvergedReason*,void*); | ||
| 32 | |||
| 33 | typedef struct { | ||
| 34 | Mat T; | ||
| 35 | Vec x1,x2,y1,y2; | ||
| 36 | PetscScalar alpha,beta,tau1,tau2,sigma; | ||
| 37 | PetscInt lastnconv; /* last value of nconv; used in stopping test */ | ||
| 38 | PetscInt nreps; /* number of repetitions of nconv; used in stopping test */ | ||
| 39 | } CTX_BRUSSEL; | ||
| 40 | |||
| 41 | 10 | int main(int argc,char **argv) | |
| 42 | { | ||
| 43 | 10 | Mat A; /* eigenvalue problem matrix */ | |
| 44 | 10 | EPS eps; /* eigenproblem solver context */ | |
| 45 | 10 | RG rg; /* region object */ | |
| 46 | 10 | PetscScalar delta1,delta2,L,h; | |
| 47 | 10 | PetscInt N=30,n,i,Istart,Iend,mpd; | |
| 48 | 10 | CTX_BRUSSEL *ctx; | |
| 49 | 10 | PetscBool terse; | |
| 50 | 10 | PetscViewer viewer; | |
| 51 | |||
| 52 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBeginUser; |
| 53 |
4/6✓ Branch 0 taken 2 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(SlepcInitialize(&argc,&argv,NULL,help)); |
| 54 | |||
| 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.
|
10 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&N,NULL)); |
| 56 |
4/6✓ Branch 0 taken 2 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(PetscPrintf(PETSC_COMM_WORLD,"\nBrusselator wave model, n=%" PetscInt_FMT "\n\n",N)); |
| 57 | |||
| 58 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 59 | Generate the matrix | ||
| 60 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 61 | |||
| 62 | /* | ||
| 63 | Create shell matrix context and set default parameters | ||
| 64 | */ | ||
| 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.
|
10 | PetscCall(PetscNew(&ctx)); |
| 66 | 10 | ctx->alpha = 2.0; | |
| 67 | 10 | ctx->beta = 5.45; | |
| 68 | 10 | delta1 = 0.008; | |
| 69 | 10 | delta2 = 0.004; | |
| 70 | 10 | L = 0.51302; | |
| 71 | |||
| 72 | /* | ||
| 73 | Look the command line for user-provided parameters | ||
| 74 | */ | ||
| 75 |
4/6✓ Branch 0 taken 2 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(PetscOptionsGetScalar(NULL,NULL,"-L",&L,NULL)); |
| 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.
|
10 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-alpha",&ctx->alpha,NULL)); |
| 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.
|
10 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-beta",&ctx->beta,NULL)); |
| 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.
|
10 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-delta1",&delta1,NULL)); |
| 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.
|
10 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-delta2",&delta2,NULL)); |
| 80 | |||
| 81 | /* | ||
| 82 | Create matrix T | ||
| 83 | */ | ||
| 84 |
4/6✓ Branch 0 taken 2 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(MatCreate(PETSC_COMM_WORLD,&ctx->T)); |
| 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.
|
10 | PetscCall(MatSetSizes(ctx->T,PETSC_DECIDE,PETSC_DECIDE,N,N)); |
| 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.
|
10 | PetscCall(MatSetFromOptions(ctx->T)); |
| 87 | |||
| 88 |
4/6✓ Branch 0 taken 2 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(MatGetOwnershipRange(ctx->T,&Istart,&Iend)); |
| 89 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1010 | for (i=Istart;i<Iend;i++) { |
| 90 |
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.
|
1000 | if (i>0) PetscCall(MatSetValue(ctx->T,i,i-1,1.0,INSERT_VALUES)); |
| 91 |
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.
|
1000 | if (i<N-1) PetscCall(MatSetValue(ctx->T,i,i+1,1.0,INSERT_VALUES)); |
| 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.
|
1000 | PetscCall(MatSetValue(ctx->T,i,i,-2.0,INSERT_VALUES)); |
| 93 | } | ||
| 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.
|
10 | PetscCall(MatAssemblyBegin(ctx->T,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.
|
10 | PetscCall(MatAssemblyEnd(ctx->T,MAT_FINAL_ASSEMBLY)); |
| 96 |
4/6✓ Branch 0 taken 2 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(MatGetLocalSize(ctx->T,&n,NULL)); |
| 97 | |||
| 98 | /* | ||
| 99 | Fill the remaining information in the shell matrix context | ||
| 100 | and create auxiliary vectors | ||
| 101 | */ | ||
| 102 | 10 | h = 1.0 / (PetscReal)(N+1); | |
| 103 | 10 | ctx->tau1 = delta1 / ((h*L)*(h*L)); | |
| 104 | 10 | ctx->tau2 = delta2 / ((h*L)*(h*L)); | |
| 105 | 10 | ctx->sigma = 0.0; | |
| 106 |
4/6✓ Branch 0 taken 2 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(VecCreateMPIWithArray(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,NULL,&ctx->x1)); |
| 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.
|
10 | PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,NULL,&ctx->x2)); |
| 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.
|
10 | PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,NULL,&ctx->y1)); |
| 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.
|
10 | PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,NULL,&ctx->y2)); |
| 110 | |||
| 111 | /* | ||
| 112 | Create the shell matrix | ||
| 113 | */ | ||
| 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.
|
10 | PetscCall(MatCreateShell(PETSC_COMM_WORLD,2*n,2*n,2*N,2*N,(void*)ctx,&A)); |
| 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.
|
10 | PetscCall(MatShellSetOperation(A,MATOP_MULT,(PetscErrorCodeFn*)MatMult_Brussel)); |
| 116 | |||
| 117 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 118 | Create the eigensolver and configure the region | ||
| 119 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 120 | |||
| 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.
|
10 | PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps)); |
| 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.
|
10 | PetscCall(EPSSetOperators(eps,A,NULL)); |
| 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.
|
10 | PetscCall(EPSSetProblemType(eps,EPS_NHEP)); |
| 124 | |||
| 125 | /* | ||
| 126 | Define the region containing the eigenvalues of interest | ||
| 127 | */ | ||
| 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.
|
10 | PetscCall(EPSGetRG(eps,&rg)); |
| 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.
|
10 | PetscCall(RGSetType(rg,RGINTERVAL)); |
| 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.
|
10 | PetscCall(RGIntervalSetEndpoints(rg,-PETSC_INFINITY,PETSC_INFINITY,-0.01,0.01)); |
| 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.
|
10 | PetscCall(RGSetComplement(rg,PETSC_TRUE)); |
| 132 | /* sort eigenvalue approximations wrt a target, otherwise convergence will be erratic */ | ||
| 133 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10 | PetscCall(EPSSetTarget(eps,0.0)); |
| 134 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10 | PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE)); |
| 135 | |||
| 136 | /* | ||
| 137 | Set solver options. In particular, we must allocate sufficient | ||
| 138 | storage for all eigenpairs that may converge (ncv). This is | ||
| 139 | application-dependent. | ||
| 140 | */ | ||
| 141 | 10 | mpd = 40; | |
| 142 |
4/6✓ Branch 0 taken 2 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(EPSSetDimensions(eps,2*mpd,3*mpd,mpd)); |
| 143 |
4/6✓ Branch 0 taken 2 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(EPSSetTolerances(eps,1e-7,2000)); |
| 144 | 10 | ctx->lastnconv = 0; | |
| 145 | 10 | ctx->nreps = 0; | |
| 146 |
4/6✓ Branch 0 taken 2 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(EPSSetStoppingTestFunction(eps,MyStoppingTest,(void*)ctx,NULL)); |
| 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.
|
10 | PetscCall(EPSSetFromOptions(eps)); |
| 148 | |||
| 149 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 150 | Solve the eigensystem and display solution | ||
| 151 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 152 | |||
| 153 |
4/6✓ Branch 0 taken 2 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(EPSSolve(eps)); |
| 154 | |||
| 155 | /* show detailed info unless -terse option is given by user */ | ||
| 156 |
4/6✓ Branch 0 taken 2 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(PetscViewerASCIIGetStdout(PETSC_COMM_WORLD,&viewer)); |
| 157 |
4/6✓ Branch 0 taken 2 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(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO_DETAIL)); |
| 158 |
4/6✓ Branch 0 taken 2 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(EPSConvergedReasonView(eps,viewer)); |
| 159 |
4/6✓ Branch 0 taken 2 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(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
| 160 |
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.
|
10 | if (!terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer)); |
| 161 |
4/6✓ Branch 0 taken 2 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(PetscViewerPopFormat(viewer)); |
| 162 | |||
| 163 |
4/6✓ Branch 0 taken 2 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(EPSDestroy(&eps)); |
| 164 |
4/6✓ Branch 0 taken 2 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(MatDestroy(&A)); |
| 165 |
4/6✓ Branch 0 taken 2 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(MatDestroy(&ctx->T)); |
| 166 |
4/6✓ Branch 0 taken 2 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(VecDestroy(&ctx->x1)); |
| 167 |
4/6✓ Branch 0 taken 2 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(VecDestroy(&ctx->x2)); |
| 168 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10 | PetscCall(VecDestroy(&ctx->y1)); |
| 169 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
10 | PetscCall(VecDestroy(&ctx->y2)); |
| 170 |
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.
|
10 | PetscCall(PetscFree(ctx)); |
| 171 |
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.
|
10 | PetscCall(SlepcFinalize()); |
| 172 | return 0; | ||
| 173 | } | ||
| 174 | |||
| 175 | 5500 | PetscErrorCode MatMult_Brussel(Mat A,Vec x,Vec y) | |
| 176 | { | ||
| 177 | 5500 | PetscInt n; | |
| 178 | 5500 | const PetscScalar *px; | |
| 179 | 5500 | PetscScalar *py; | |
| 180 | 5500 | CTX_BRUSSEL *ctx; | |
| 181 | |||
| 182 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
5500 | PetscFunctionBeginUser; |
| 183 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(MatShellGetContext(A,&ctx)); |
| 184 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(MatGetLocalSize(ctx->T,&n,NULL)); |
| 185 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecGetArrayRead(x,&px)); |
| 186 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecGetArray(y,&py)); |
| 187 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecPlaceArray(ctx->x1,px)); |
| 188 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecPlaceArray(ctx->x2,px+n)); |
| 189 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecPlaceArray(ctx->y1,py)); |
| 190 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecPlaceArray(ctx->y2,py+n)); |
| 191 | |||
| 192 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(MatMult(ctx->T,ctx->x1,ctx->y1)); |
| 193 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecScale(ctx->y1,ctx->tau1)); |
| 194 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecAXPY(ctx->y1,ctx->beta - 1.0 + ctx->sigma,ctx->x1)); |
| 195 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecAXPY(ctx->y1,ctx->alpha * ctx->alpha,ctx->x2)); |
| 196 | |||
| 197 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(MatMult(ctx->T,ctx->x2,ctx->y2)); |
| 198 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecScale(ctx->y2,ctx->tau2)); |
| 199 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecAXPY(ctx->y2,-ctx->beta,ctx->x1)); |
| 200 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecAXPY(ctx->y2,-ctx->alpha * ctx->alpha + ctx->sigma,ctx->x2)); |
| 201 | |||
| 202 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecRestoreArrayRead(x,&px)); |
| 203 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecRestoreArray(y,&py)); |
| 204 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecResetArray(ctx->x1)); |
| 205 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecResetArray(ctx->x2)); |
| 206 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecResetArray(ctx->y1)); |
| 207 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5500 | PetscCall(VecResetArray(ctx->y2)); |
| 208 |
5/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 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
1100 | PetscFunctionReturn(PETSC_SUCCESS); |
| 209 | } | ||
| 210 | |||
| 211 | /* | ||
| 212 | Function for user-defined stopping test. | ||
| 213 | |||
| 214 | Ignores the value of nev. It only takes into account the number of | ||
| 215 | eigenpairs that have converged in recent outer iterations (restarts); | ||
| 216 | if no new eigenvalues have converged in the last few restarts, | ||
| 217 | we stop the iteration, assuming that no more eigenvalues are present | ||
| 218 | inside the region. | ||
| 219 | */ | ||
| 220 | 260 | PetscErrorCode MyStoppingTest(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ptr) | |
| 221 | { | ||
| 222 | 260 | CTX_BRUSSEL *ctx = (CTX_BRUSSEL*)ptr; | |
| 223 | |||
| 224 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
260 | PetscFunctionBeginUser; |
| 225 | /* check usual termination conditions, but ignoring the case nconv>=nev */ | ||
| 226 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
260 | PetscCall(EPSStoppingBasic(eps,its,max_it,nconv,PETSC_INT_MAX,reason,NULL)); |
| 227 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
260 | if (*reason==EPS_CONVERGED_ITERATING) { |
| 228 | /* check if nconv is the same as before */ | ||
| 229 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
260 | if (nconv==ctx->lastnconv) ctx->nreps++; |
| 230 | else { | ||
| 231 | 45 | ctx->lastnconv = nconv; | |
| 232 | 45 | ctx->nreps = 0; | |
| 233 | } | ||
| 234 | /* check if no eigenvalues converged in last 10 restarts */ | ||
| 235 |
4/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 10 times.
|
260 | if (nconv && ctx->nreps>10) *reason = EPS_CONVERGED_USER; |
| 236 | } | ||
| 237 |
5/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 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
52 | PetscFunctionReturn(PETSC_SUCCESS); |
| 238 | } | ||
| 239 | |||
| 240 | /*TEST | ||
| 241 | |||
| 242 | test: | ||
| 243 | suffix: 1 | ||
| 244 | args: -n 100 -terse | ||
| 245 | requires: !single | ||
| 246 | |||
| 247 | TEST*/ | ||
| 248 |