Actual source code: filter.h

  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: typedef struct {
 14:   PetscErrorCode  (*computeoperator)(ST,Mat*);
 15:   PetscErrorCode  (*getthreshold)(ST,PetscReal*);
 16:   PetscErrorCode  (*destroy)(ST);

 18:   STFilterType    type;            /* the filter method */
 19:   PetscReal       inta,intb;       /* bounds of the interval of desired eigenvalues */
 20:   PetscReal       left,right;      /* approximate left and right bounds of the interval containing all eigenvalues */
 21:   PetscInt        polyDegree;      /* degree of the polynomial filter */
 22:   STFilterDamping damping;         /* the type of damping */
 23:   Mat             T;               /* the matrix used to build the filter */
 24:   PetscBool       filtch;          /* filter parameters have changed since last setup */
 25:   Mat             *W;              /* work matrices for the matrix-matrix application of the filter */
 26:   PetscInt        nW;              /* number of W matrices */
 27:   void            *data;           /* placeholder for method-specific stuff */
 28: } ST_FILTER;

 30: /* FILTLAN */

 32: /* IntervalOptions structure used by GetIntervals */
 33: struct _n_FILTLAN_IOP {
 34:   PetscReal   intervalWeights[5];  /* weight of the subintervals (5 in mid-pass, 3 in high-pass) */
 35:   PetscReal   transIntervalRatio;  /* the (relative) length of transition interval */
 36:   PetscBool   reverseInterval;     /* whether to reverse the interval or not; effective only for mid-pass filters */
 37:   PetscReal   initialPlateau;      /* initial plateau relative to the length of interval */
 38:   PetscReal   plateauShrinkRate;   /* the rate at which the plateau shrinks at each iteration */
 39:   PetscReal   initialShiftStep;    /* initial shift step relative to the length of interval */
 40:   PetscReal   shiftStepExpanRate;  /* the rate at which the shift step expands */
 41:   PetscInt    maxInnerIter;        /* maximum number of inner iterations to determine the (transition) intervals */
 42:   PetscReal   yLimitTol;           /* tolerance allowed for |p(inta)-p(intb)| in a mid-pass filter p(x) */
 43:   PetscInt    maxOuterIter;        /* maximum number of outer iterations to determine the (transition) intervals */
 44:   PetscInt    numGridPoints;       /* number of grid points, used to measure the maximum p(z) for z not in the interval*/
 45:   PetscReal   yBottomLine;         /* p(x) should be greater than this value for x in the interval */
 46:   PetscReal   yRippleLimit;        /* limit of height of ripples relative to the bottom of polynomial values */
 47: };
 48: typedef struct _n_FILTLAN_IOP *FILTLAN_IOP;

 50: /* PolynomialFilterInfo structure filled by GetIntervals */
 51: struct _n_FILTLAN_PFI {
 52:   PetscInt    filterType;          /* 1 = high-pass filter (or low-pass filter with conversion); 2 = mid-pass filter */
 53:   PetscInt    filterOK;            /* 0 = no acceptable found; 1 = OK filter found; 2 = optimal filter found */
 54:   PetscReal   filterQualityIndex;  /* between 0.0 and 1.0; the higher the better quality of the filter */
 55:   PetscInt    numIter;             /* number of iterations to get the (transition) intervals */
 56:   PetscInt    totalNumIter;        /* total number of iterations performed */
 57:   PetscReal   yLimit;              /* lowest polynomial value P(z) for z in the interval [a1,b1] of desired eigenvalues */
 58:   PetscReal   ySummit;             /* height of (highest, if more than one) summit in interval [a1,b1] of desired evals */
 59:   PetscInt    numLeftSteps;        /* number of steps moving leftward */
 60:   PetscReal   yLeftSummit;         /* height of highest summit in the left-hand side of the interval of desired evals */
 61:   PetscReal   yLeftBottom;         /* height of lowest bottom in the left-hand side of the interval of desired evals */
 62:   PetscReal   yLimitGap;           /* |p(a1)-p(b1)|, where [a1,b1] is the interval of desired eigenvalues */
 63:   PetscInt    numRightSteps;       /* number of steps moving rightward */
 64:   PetscReal   yRightSummit;        /* height of highest summit in the right-hand side of the interval of desired evals */
 65:   PetscReal   yRightBottom;        /* height of lowest bottom in the right-hand side of the interval of desired evals */
 66: };
 67: typedef struct _n_FILTLAN_PFI *FILTLAN_PFI;

 69: struct _n_FILTLAN_CTX {
 70:   PetscInt      baseDegree;        /* left and right degrees of the base filter for each interval */
 71:   FILTLAN_IOP   opts;              /* interval options */
 72:   FILTLAN_PFI   filterInfo;        /* polynomial filter info */
 73:   PetscReal     frame[4];          /* outer and inner intervals:
 74:                                       [frame[0],frame[3]] (tightly) contains all eigenvalues
 75:                                       [frame[1],frame[2]] is the interval in which the eigenvalues are sought */
 76:   PetscReal     intervals[6];      /* the range of eigenvalues is partitioned into intervals which determine
 77:                                       the base filter approximated by a polynomial filter;
 78:                                       the j-th interval is [intervals(j),intervals(j+1)) */
 79:   PetscReal     intervals2[6];     /* modified intervals */
 80:   PetscReal     *baseFilter;       /* coefficients of the base filter (a piecewise polynomial) */
 81: };
 82: typedef struct _n_FILTLAN_CTX *FILTLAN_CTX;

 84: SLEPC_INTERN PetscErrorCode STCreate_Filter_FILTLAN(ST);

 86: /* Chebyshev series */

 88: struct _n_CHEBYSHEV_CTX {
 89:   PetscReal   *damping_coeffs;
 90:   PetscReal   *coeffs;
 91: };
 92: typedef struct _n_CHEBYSHEV_CTX *CHEBYSHEV_CTX;

 94: SLEPC_INTERN PetscErrorCode STCreate_Filter_Chebyshev(ST);