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