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 | This example implements one of the problems found at | ||
12 | NLEVP: A Collection of Nonlinear Eigenvalue Problems, | ||
13 | The University of Manchester. | ||
14 | The details of the collection can be found at: | ||
15 | [1] T. Betcke et al., "NLEVP: A Collection of Nonlinear Eigenvalue | ||
16 | Problems", ACM Trans. Math. Software 39(2), Article 7, 2013. | ||
17 | |||
18 | The pdde_stability problem is a complex-symmetric QEP from the stability | ||
19 | analysis of a discretized partial delay-differential equation. It requires | ||
20 | complex scalars. | ||
21 | */ | ||
22 | |||
23 | static char help[] = "Stability analysis of a discretized partial delay-differential equation.\n\n" | ||
24 | "The command line options are:\n" | ||
25 | " -m <m>, grid size, the matrices have dimension n=m*m.\n" | ||
26 | " -c <a0,b0,a1,b1,a2,b2,phi1>, comma-separated list of 7 real parameters.\n\n"; | ||
27 | |||
28 | #include <slepcpep.h> | ||
29 | |||
30 | #define NMAT 3 | ||
31 | |||
32 | /* | ||
33 | Function for user-defined eigenvalue ordering criterion. | ||
34 | |||
35 | Given two eigenvalues ar+i*ai and br+i*bi, the subroutine must choose | ||
36 | one of them as the preferred one according to the criterion. | ||
37 | In this example, the preferred value is the one with absolute value closest to 1. | ||
38 | */ | ||
39 | 83150 | PetscErrorCode MyEigenSort(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *r,void *ctx) | |
40 | { | ||
41 | 83150 | PetscReal aa,ab; | |
42 | |||
43 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
83150 | PetscFunctionBeginUser; |
44 | 83150 | aa = PetscAbsReal(SlepcAbsEigenvalue(ar,ai)-PetscRealConstant(1.0)); | |
45 | 83150 | ab = PetscAbsReal(SlepcAbsEigenvalue(br,bi)-PetscRealConstant(1.0)); | |
46 |
3/4✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 5 times.
|
83150 | *r = aa > ab ? 1 : (aa < ab ? -1 : 0); |
47 |
5/12✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 1 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 1 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
83150 | PetscFunctionReturn(PETSC_SUCCESS); |
48 | } | ||
49 | |||
50 | 15 | int main(int argc,char **argv) | |
51 | { | ||
52 | 15 | Mat A[NMAT]; /* problem matrices */ | |
53 | 15 | PEP pep; /* polynomial eigenproblem solver context */ | |
54 | 15 | PetscInt m=15,n,II,Istart,Iend,i,j,k; | |
55 | 15 | PetscReal h,xi,xj,c[7] = { 2, .3, -2, .2, -2, -.3, -PETSC_PI/2 }; | |
56 | 15 | PetscScalar alpha,beta,gamma; | |
57 | 15 | PetscBool flg,terse; | |
58 | |||
59 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
15 | PetscFunctionBeginUser; |
60 |
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.
|
15 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
61 | #if !defined(PETSC_USE_COMPLEX) | ||
62 | SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"This example requires complex scalars"); | ||
63 | #endif | ||
64 | |||
65 |
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.
|
15 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL)); |
66 | 15 | n = m*m; | |
67 | 15 | h = PETSC_PI/(m+1); | |
68 | 15 | gamma = PetscExpScalar(PETSC_i*c[6]); | |
69 | 15 | gamma = gamma/PetscAbsScalar(gamma); | |
70 | 15 | k = 7; | |
71 |
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.
|
15 | PetscCall(PetscOptionsGetRealArray(NULL,NULL,"-c",c,&k,&flg)); |
72 |
1/6✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
15 | PetscCheck(!flg || k==7,PETSC_COMM_WORLD,PETSC_ERR_USER,"The number of parameters -c should be 7, you provided %" PetscInt_FMT,k); |
73 |
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.
|
15 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nPDDE stability, n=%" PetscInt_FMT " (m=%" PetscInt_FMT ")\n\n",n,m)); |
74 | |||
75 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
76 | Compute the polynomial matrices | ||
77 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
78 | |||
79 | /* initialize matrices */ | ||
80 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
60 | for (i=0;i<NMAT;i++) { |
81 |
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.
|
45 | PetscCall(MatCreate(PETSC_COMM_WORLD,&A[i])); |
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.
|
45 | PetscCall(MatSetSizes(A[i],PETSC_DECIDE,PETSC_DECIDE,n,n)); |
83 |
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.
|
45 | PetscCall(MatSetFromOptions(A[i])); |
84 | } | ||
85 |
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.
|
15 | PetscCall(MatGetOwnershipRange(A[0],&Istart,&Iend)); |
86 | |||
87 | /* A[1] has a pattern similar to the 2D Laplacian */ | ||
88 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
3390 | for (II=Istart;II<Iend;II++) { |
89 | 3375 | i = II/m; j = II-i*m; | |
90 | 3375 | xi = (i+1)*h; xj = (j+1)*h; | |
91 | 3375 | alpha = c[0]+c[1]*PetscSinReal(xi)+gamma*(c[2]+c[3]*xi*(1.0-PetscExpReal(xi-PETSC_PI))); | |
92 | 3375 | beta = c[0]+c[1]*PetscSinReal(xj)-gamma*(c[2]+c[3]*xj*(1.0-PetscExpReal(xj-PETSC_PI))); | |
93 |
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.
|
3375 | PetscCall(MatSetValue(A[1],II,II,alpha+beta-4.0/(h*h),INSERT_VALUES)); |
94 |
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.
|
3375 | if (j>0) PetscCall(MatSetValue(A[1],II,II-1,1.0/(h*h),INSERT_VALUES)); |
95 |
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.
|
3375 | if (j<m-1) PetscCall(MatSetValue(A[1],II,II+1,1.0/(h*h),INSERT_VALUES)); |
96 |
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.
|
3375 | if (i>0) PetscCall(MatSetValue(A[1],II,II-m,1.0/(h*h),INSERT_VALUES)); |
97 |
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.
|
3375 | if (i<m-1) PetscCall(MatSetValue(A[1],II,II+m,1.0/(h*h),INSERT_VALUES)); |
98 | } | ||
99 | |||
100 | /* A[0] and A[2] are diagonal */ | ||
101 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
3390 | for (II=Istart;II<Iend;II++) { |
102 | 3375 | i = II/m; j = II-i*m; | |
103 | 3375 | xi = (i+1)*h; xj = (j+1)*h; | |
104 | 3375 | alpha = c[4]+c[5]*xi*(PETSC_PI-xi); | |
105 | 3375 | beta = c[4]+c[5]*xj*(PETSC_PI-xj); | |
106 |
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.
|
3375 | PetscCall(MatSetValue(A[0],II,II,alpha,INSERT_VALUES)); |
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.
|
3375 | PetscCall(MatSetValue(A[2],II,II,beta,INSERT_VALUES)); |
108 | } | ||
109 | |||
110 | /* assemble matrices */ | ||
111 |
7/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 5 times.
✓ Branch 3 taken 4 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 1 times.
|
60 | for (i=0;i<NMAT;i++) PetscCall(MatAssemblyBegin(A[i],MAT_FINAL_ASSEMBLY)); |
112 |
7/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 5 times.
✓ Branch 3 taken 4 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 1 times.
|
60 | for (i=0;i<NMAT;i++) PetscCall(MatAssemblyEnd(A[i],MAT_FINAL_ASSEMBLY)); |
113 | |||
114 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
115 | Create the eigensolver and solve the problem | ||
116 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
117 | |||
118 |
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.
|
15 | PetscCall(PEPCreate(PETSC_COMM_WORLD,&pep)); |
119 |
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.
|
15 | PetscCall(PEPSetOperators(pep,NMAT,A)); |
120 |
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.
|
15 | PetscCall(PEPSetEigenvalueComparison(pep,MyEigenSort,NULL)); |
121 |
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.
|
15 | PetscCall(PEPSetDimensions(pep,4,PETSC_DETERMINE,PETSC_DETERMINE)); |
122 |
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.
|
15 | PetscCall(PEPSetFromOptions(pep)); |
123 |
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.
|
15 | PetscCall(PEPSolve(pep)); |
124 | |||
125 | /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
126 | Display solution and clean up | ||
127 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ | ||
128 | |||
129 | /* show detailed info unless -terse option is given by user */ | ||
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.
|
15 | PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse)); |
131 |
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.
|
15 | if (terse) PetscCall(PEPErrorView(pep,PEP_ERROR_BACKWARD,NULL)); |
132 | else { | ||
133 | ✗ | PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL)); | |
134 | ✗ | PetscCall(PEPConvergedReasonView(pep,PETSC_VIEWER_STDOUT_WORLD)); | |
135 | ✗ | PetscCall(PEPErrorView(pep,PEP_ERROR_BACKWARD,PETSC_VIEWER_STDOUT_WORLD)); | |
136 | ✗ | PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); | |
137 | } | ||
138 |
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.
|
15 | PetscCall(PEPDestroy(&pep)); |
139 |
7/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 5 times.
✓ Branch 3 taken 4 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 1 times.
|
60 | for (i=0;i<NMAT;i++) PetscCall(MatDestroy(&A[i])); |
140 |
3/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
15 | PetscCall(SlepcFinalize()); |
141 | return 0; | ||
142 | } | ||
143 | |||
144 | /*TEST | ||
145 | |||
146 | build: | ||
147 | requires: complex | ||
148 | |||
149 | test: | ||
150 | suffix: 1 | ||
151 | args: -pep_type {{toar qarnoldi linear}} -pep_ncv 25 -terse | ||
152 | requires: complex double | ||
153 | |||
154 | TEST*/ | ||
155 |