| 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[] = "Use the matrix exponential to compute rightmost eigenvalues.\n\n" | ||
| 12 | "Same problem as ex9.c but with explicitly created matrix. The command line options are:\n" | ||
| 13 | " -n <n>, where <n> = block dimension of the 2x2 block matrix.\n" | ||
| 14 | " -L <L>, where <L> = bifurcation parameter.\n" | ||
| 15 | " -alpha <alpha>, -beta <beta>, -delta1 <delta1>, -delta2 <delta2>,\n" | ||
| 16 | " where <alpha> <beta> <delta1> <delta2> = model parameters.\n\n"; | ||
| 17 | |||
| 18 | #include <slepceps.h> | ||
| 19 | #include <slepcmfn.h> | ||
| 20 | |||
| 21 | /* | ||
| 22 | This example computes the eigenvalues with largest real part of the | ||
| 23 | following matrix | ||
| 24 | |||
| 25 | A = [ tau1*T+(beta-1)*I alpha^2*I | ||
| 26 | -beta*I tau2*T-alpha^2*I ], | ||
| 27 | |||
| 28 | where | ||
| 29 | |||
| 30 | T = tridiag{1,-2,1} | ||
| 31 | h = 1/(n+1) | ||
| 32 | tau1 = delta1/(h*L)^2 | ||
| 33 | tau2 = delta2/(h*L)^2 | ||
| 34 | |||
| 35 | but it builds A explicitly, as opposed to ex9.c | ||
| 36 | */ | ||
| 37 | |||
| 38 | /* Routines for shell spectral transformation */ | ||
| 39 | PetscErrorCode STApply_Exp(ST,Vec,Vec); | ||
| 40 | PetscErrorCode STBackTransform_Exp(ST,PetscInt,PetscScalar*,PetscScalar*); | ||
| 41 | |||
| 42 | 60 | int main(int argc,char **argv) | |
| 43 | { | ||
| 44 | 60 | Mat A; /* operator matrix */ | |
| 45 | 60 | EPS eps; /* eigenproblem solver context */ | |
| 46 | 60 | ST st; /* spectral transformation context */ | |
| 47 | 60 | MFN mfn; /* matrix function solver object to compute exp(A)*v */ | |
| 48 | 60 | FN f; | |
| 49 | 60 | EPSType type; | |
| 50 | 60 | PetscScalar alpha,beta,tau1,tau2,delta1,delta2,L,h; | |
| 51 | 60 | PetscInt n=30,i,Istart,Iend,nev; | |
| 52 | 60 | PetscBool isShell,terse; | |
| 53 | |||
| 54 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
60 | PetscFunctionBeginUser; |
| 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.
|
60 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
| 56 | #if defined(PETSC_HAVE_COMPLEX) | ||
| 57 |
4/6✓ Branch 0 taken 2 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(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); |
| 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(PetscPrintf(PETSC_COMM_WORLD,"\nBrusselator wave model with matrix exponential, n=%" PetscInt_FMT "\n\n",n)); |
| 59 | |||
| 60 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 61 | Generate the matrix | ||
| 62 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 63 | |||
| 64 | 60 | alpha = 2.0; | |
| 65 | 60 | beta = 5.45; | |
| 66 | 60 | delta1 = 0.008; | |
| 67 | 60 | delta2 = 0.004; | |
| 68 | 60 | L = 0.51302; | |
| 69 | |||
| 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.
|
60 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-L",&L,NULL)); |
| 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.
|
60 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-alpha",&alpha,NULL)); |
| 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.
|
60 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-beta",&beta,NULL)); |
| 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.
|
60 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-delta1",&delta1,NULL)); |
| 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.
|
60 | PetscCall(PetscOptionsGetScalar(NULL,NULL,"-delta2",&delta2,NULL)); |
| 75 | |||
| 76 | 60 | h = 1.0 / (PetscReal)(n+1); | |
| 77 | 60 | tau1 = delta1 / ((h*L)*(h*L)); | |
| 78 | 60 | tau2 = delta2 / ((h*L)*(h*L)); | |
| 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.
|
60 | PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); |
| 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.
|
60 | PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,2*n,2*n)); |
| 82 |
4/6✓ Branch 0 taken 2 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(MatSetFromOptions(A)); |
| 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.
|
60 | PetscCall(MatGetOwnershipRange(A,&Istart,&Iend)); |
| 85 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
7140 | for (i=Istart;i<Iend;i++) { |
| 86 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
7080 | if (i<n) { /* upper blocks */ |
| 87 |
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.
|
3540 | if (i>0) PetscCall(MatSetValue(A,i,i-1,tau1,INSERT_VALUES)); |
| 88 |
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.
|
3540 | if (i<n-1) PetscCall(MatSetValue(A,i,i+1,tau1,INSERT_VALUES)); |
| 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.
|
3540 | PetscCall(MatSetValue(A,i,i,-2.0*tau1+beta-1.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.
|
3540 | PetscCall(MatSetValue(A,i,i+n,alpha*alpha,INSERT_VALUES)); |
| 91 | } else { /* lower blocks */ | ||
| 92 |
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.
|
3540 | if (i>n) PetscCall(MatSetValue(A,i,i-1,tau2,INSERT_VALUES)); |
| 93 |
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.
|
3540 | if (i<2*n-1) PetscCall(MatSetValue(A,i,i+1,tau2,INSERT_VALUES)); |
| 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.
|
3540 | PetscCall(MatSetValue(A,i,i,-2.0*tau2-alpha*alpha,INSERT_VALUES)); |
| 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.
|
7080 | PetscCall(MatSetValue(A,i,i-n,-beta,INSERT_VALUES)); |
| 96 | } | ||
| 97 | } | ||
| 98 |
4/6✓ Branch 0 taken 2 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(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); |
| 99 |
4/6✓ Branch 0 taken 2 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(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); |
| 100 | |||
| 101 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 102 | Create the eigensolver and set various options | ||
| 103 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 104 | |||
| 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.
|
60 | PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps)); |
| 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.
|
60 | PetscCall(EPSSetOperators(eps,A,NULL)); |
| 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.
|
60 | PetscCall(EPSSetProblemType(eps,EPS_NHEP)); |
| 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.
|
60 | PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL)); |
| 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.
|
60 | PetscCall(EPSGetST(eps,&st)); |
| 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.
|
60 | PetscCall(STSetType(st,STSHELL)); |
| 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.
|
60 | PetscCall(EPSSetFromOptions(eps)); |
| 112 | |||
| 113 | /* | ||
| 114 | Initialize shell spectral transformation | ||
| 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.
|
60 | PetscCall(PetscObjectTypeCompare((PetscObject)st,STSHELL,&isShell)); |
| 117 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
60 | if (isShell) { |
| 118 | |||
| 119 | /* Create the MFN object to be used by the spectral transform */ | ||
| 120 |
4/6✓ Branch 0 taken 2 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(MFNCreate(PETSC_COMM_WORLD,&mfn)); |
| 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(MFNSetOperator(mfn,A)); |
| 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(MFNGetFN(mfn,&f)); |
| 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(FNSetType(f,FNEXP)); |
| 124 |
4/6✓ Branch 0 taken 2 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(FNSetScale(f,0.03,1.0)); /* this can be set with -fn_scale */ |
| 125 |
4/6✓ Branch 0 taken 2 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(MFNSetFromOptions(mfn)); |
| 126 | |||
| 127 | /* Set callback functions */ | ||
| 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(STShellSetApply(st,STApply_Exp)); |
| 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(STShellSetBackTransform(st,STBackTransform_Exp)); |
| 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(STShellSetContext(st,mfn)); |
| 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(PetscObjectSetName((PetscObject)st,"STEXP")); |
| 132 | } | ||
| 133 | |||
| 134 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 135 | Solve the eigensystem | ||
| 136 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 137 | |||
| 138 |
4/6✓ Branch 0 taken 2 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(EPSSolve(eps)); |
| 139 |
4/6✓ Branch 0 taken 2 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(EPSGetType(eps,&type)); |
| 140 |
4/6✓ Branch 0 taken 2 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(PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type)); |
| 141 |
4/6✓ Branch 0 taken 2 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(EPSGetDimensions(eps,&nev,NULL,NULL)); |
| 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.
|
60 | PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of requested eigenvalues: %" PetscInt_FMT "\n",nev)); |
| 143 | |||
| 144 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 145 | Display solution and clean up | ||
| 146 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
| 147 | |||
| 148 | /* show detailed info unless -terse option is given by user */ | ||
| 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.
|
60 | PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
| 150 |
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.
|
60 | if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL)); |
| 151 | else { | ||
| 152 | ✗ | PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL)); | |
| 153 | ✗ | PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD)); | |
| 154 | ✗ | PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD)); | |
| 155 | ✗ | PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); | |
| 156 | } | ||
| 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.
|
60 | PetscCall(EPSDestroy(&eps)); |
| 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.
|
60 | PetscCall(MatDestroy(&A)); |
| 159 |
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.
|
60 | if (isShell) PetscCall(MFNDestroy(&mfn)); |
| 160 | #else | ||
| 161 | SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"This example requires C99 complex numbers"); | ||
| 162 | #endif | ||
| 163 |
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.
|
60 | PetscCall(SlepcFinalize()); |
| 164 | return 0; | ||
| 165 | } | ||
| 166 | |||
| 167 | /* ------------------------------------------------------------------- */ | ||
| 168 | /* | ||
| 169 | STBackTransform_Exp - Undoes the exp(A) transformation by taking logarithms. | ||
| 170 | |||
| 171 | Input Parameters: | ||
| 172 | + st - spectral transformation context | ||
| 173 | - n - number of eigenvalues to transform | ||
| 174 | |||
| 175 | Input/Output Parameters: | ||
| 176 | + eigr - pointer to real part of eigenvalues | ||
| 177 | - eigi - pointer to imaginary part of eigenvalues | ||
| 178 | */ | ||
| 179 | 5890 | PetscErrorCode STBackTransform_Exp(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi) | |
| 180 | { | ||
| 181 | #if defined(PETSC_HAVE_COMPLEX) | ||
| 182 | 5890 | PetscInt j; | |
| 183 | 5890 | MFN mfn; | |
| 184 | 5890 | FN fn; | |
| 185 | 5890 | PetscScalar tau,eta; | |
| 186 | #if !defined(PETSC_USE_COMPLEX) | ||
| 187 | 1785 | PetscComplex theta,lambda; | |
| 188 | #endif | ||
| 189 | |||
| 190 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
5890 | PetscFunctionBeginUser; |
| 191 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
5890 | PetscCall(STShellGetContext(st,&mfn)); |
| 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.
|
5890 | PetscCall(MFNGetFN(mfn,&fn)); |
| 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.
|
5890 | PetscCall(FNGetScale(fn,&tau,&eta)); |
| 194 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
17710 | for (j=0;j<n;j++) { |
| 195 | #if defined(PETSC_USE_COMPLEX) | ||
| 196 | 8230 | eigr[j] = PetscLogComplex(eigr[j]/eta)/tau; | |
| 197 | #else | ||
| 198 | 3590 | theta = PetscCMPLX(eigr[j],eigi[j])/eta; | |
| 199 | 3590 | lambda = PetscLogComplex(theta)/tau; | |
| 200 | 3590 | eigr[j] = PetscRealPartComplex(lambda); | |
| 201 | 3590 | eigi[j] = PetscImaginaryPartComplex(lambda); | |
| 202 | #endif | ||
| 203 | } | ||
| 204 |
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.
|
1178 | PetscFunctionReturn(PETSC_SUCCESS); |
| 205 | #else | ||
| 206 | return 0; | ||
| 207 | #endif | ||
| 208 | } | ||
| 209 | |||
| 210 | /* | ||
| 211 | STApply_Exp - Applies the operator exp(tau*A) to a given vector using an MFN object. | ||
| 212 | |||
| 213 | Input Parameters: | ||
| 214 | + st - spectral transformation context | ||
| 215 | - x - input vector | ||
| 216 | |||
| 217 | Output Parameter: | ||
| 218 | . y - output vector | ||
| 219 | */ | ||
| 220 | 570 | PetscErrorCode STApply_Exp(ST st,Vec x,Vec y) | |
| 221 | { | ||
| 222 | 570 | MFN mfn; | |
| 223 | |||
| 224 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
570 | PetscFunctionBeginUser; |
| 225 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
570 | PetscCall(STShellGetContext(st,&mfn)); |
| 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.
|
570 | PetscCall(MFNSolve(mfn,x,y)); |
| 227 |
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.
|
114 | PetscFunctionReturn(PETSC_SUCCESS); |
| 228 | } | ||
| 229 | |||
| 230 | /*TEST | ||
| 231 | |||
| 232 | testset: | ||
| 233 | args: -eps_nev 4 -mfn_ncv 16 -terse | ||
| 234 | requires: c99_complex !single | ||
| 235 | filter: sed -e "s/-2/+2/g" | ||
| 236 | output_file: output/ex36_1.out | ||
| 237 | test: | ||
| 238 | suffix: 1 | ||
| 239 | requires: !__float128 | ||
| 240 | test: | ||
| 241 | suffix: 1_quad | ||
| 242 | args: -eps_tol 1e-14 | ||
| 243 | requires: __float128 | ||
| 244 | |||
| 245 | test: | ||
| 246 | suffix: 2 | ||
| 247 | args: -n 56 -eps_nev 2 -st_type sinvert -eps_target -390 -eps_target_magnitude -eps_type power | ||
| 248 | args: -eps_power_shift_type {{constant rayleigh}} -eps_two_sided {{0 1}} -eps_tol 1e-14 -terse | ||
| 249 | requires: c99_complex !single | ||
| 250 | filter: sed -e "s/[+-]0\.0*i//g" | ||
| 251 | |||
| 252 | test: | ||
| 253 | suffix: 3 | ||
| 254 | args: -n 100 -st_type sinvert -eps_type ciss -rg_type ellipse -rg_ellipse_center 0 -rg_ellipse_radius 6 -eps_all -eps_tol 1e-6 -terse | ||
| 255 | requires: c99_complex !single | ||
| 256 | filter: sed -e "s/-3.37036-3.55528i, -3.37036+3.55528i/-3.37036+3.55528i, -3.37036-3.55528i/" -e "s/-1.79853-3.03216i, -1.79853+3.03216i/-1.79853+3.03216i, -1.79853-3.03216i/" -e "s/-0.67471-2.52856i, -0.67471+2.52856i/-0.67471+2.52856i, -0.67471-2.52856i/" -e "s/0.00002-2.13950i, 0.00002+2.13950i/0.00002+2.13950i, 0.00002-2.13950i/" | ||
| 257 | |||
| 258 | TEST*/ | ||
| 259 |