1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #if !defined(SLEPC_FILTER_H)
12: #define SLEPC_FILTER_H 14: /* IntervalOptions structure used by GetIntervals */
15: struct _n_FILTLAN_IOP {
16: PetscReal intervalWeights[5]; /* weight of the subintervals (5 in mid-pass, 3 in high-pass) */
17: PetscReal transIntervalRatio; /* the (relative) length of transition interval */
18: PetscBool reverseInterval; /* whether to reverse the interval or not; effective only for mid-pass filters */
19: PetscReal initialPlateau; /* initial plateau relative to the length of interval */
20: PetscReal plateauShrinkRate; /* the rate at which the plateau shrinks at each iteration */
21: PetscReal initialShiftStep; /* initial shift step relative to the length of interval */
22: PetscReal shiftStepExpanRate; /* the rate at which the shift step expands */
23: PetscInt maxInnerIter; /* maximum number of inner iterations to determine the (transition) intervals */
24: PetscReal yLimitTol; /* tolerance allowed for |p(inta)-p(intb)| in a mid-pass filter p(x) */
25: PetscInt maxOuterIter; /* maximum number of outer iterations to determine the (transition) intervals */
26: PetscInt numGridPoints; /* number of grid points, used to measure the maximum p(z) for z not in the interval*/
27: PetscReal yBottomLine; /* p(x) should be greater than this value for x in the interval */
28: PetscReal yRippleLimit; /* limit of height of ripples relative to the bottom of polynomial values */
29: };
30: typedef struct _n_FILTLAN_IOP *FILTLAN_IOP;
32: /* PolynomialFilterInfo structure filled by GetIntervals */
33: struct _n_FILTLAN_PFI {
34: PetscInt filterType; /* 1 = high-pass filter (or low-pass filter with conversion); 2 = mid-pass filter */
35: PetscInt filterOK; /* 0 = no acceptable found; 1 = OK filter found; 2 = optimal filter found */
36: PetscReal filterQualityIndex; /* between 0.0 and 1.0; the higher the better quality of the filter */
37: PetscInt numIter; /* number of iterations to get the (transition) intervals */
38: PetscInt totalNumIter; /* total number of iterations performed */
39: PetscReal yLimit; /* lowest polynomial value P(z) for z in the interval [a1,b1] of desired eigenvalues */
40: PetscReal ySummit; /* height of (highest, if more than one) summit in interval [a1,b1] of desired evals */
41: PetscInt numLeftSteps; /* number of steps moving leftward */
42: PetscReal yLeftSummit; /* height of highest summit in the left-hand side of the interval of desired evals */
43: PetscReal yLeftBottom; /* height of lowest bottom in the left-hand side of the interval of desired evals */
44: PetscReal yLimitGap; /* |p(a1)-p(b1)|, where [a1,b1] is the interval of desired eigenvalues */
45: PetscInt numRightSteps; /* number of steps moving rightward */
46: PetscReal yRightSummit; /* height of highest summit in the right-hand side of the interval of desired evals */
47: PetscReal yRightBottom; /* height of lowest bottom in the right-hand side of the interval of desired evals */
48: };
49: typedef struct _n_FILTLAN_PFI *FILTLAN_PFI;
51: typedef struct {
52: /* user options */
53: PetscReal inta,intb; /* bounds of the interval of desired eigenvalues */
54: PetscReal left,right; /* approximate left and right bounds of the interval containing all eigenvalues */
55: PetscInt polyDegree; /* degree of s(z), with z*s(z) the polynomial filter */
56: PetscInt baseDegree; /* left and right degrees of the base filter for each interval */
57: FILTLAN_IOP opts; /* interval options */
58: /* internal variables */
59: FILTLAN_PFI filterInfo; /* polynomial filter info */
60: PetscReal frame[4]; /* outer and inner intervals:
61: [frame[0],frame[3]] (tightly) contains all eigenvalues
62: [frame[1],frame[2]] is the interval in which the eigenvalues are sought */
63: PetscReal intervals[6]; /* the range of eigenvalues is partitioned into intervals which determine
64: the base filter approximated by a polynomial filter;
65: the j-th interval is [intervals(j),intervals(j+1)) */
66: PetscReal intervals2[6]; /* modified intervals */
67: PetscReal *baseFilter; /* coefficients of the base filter (a piecewise polynomial) */
68: } ST_FILTER;
70: SLEPC_INTERN PetscErrorCode STFilter_FILTLAN_Apply(ST,Vec,Vec);
71: SLEPC_INTERN PetscErrorCode STFilter_FILTLAN_setFilter(ST);
73: #endif