Actual source code: filter.c
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: */
10: /*
11: Filter spectral transformation, to encapsulate polynomial filters
12: */
14: #include <slepc/private/stimpl.h>
15: #include "filter.h"
17: const char *STFilterTypes[] = {"","FILTLAN","CHEBYSHEV","STFilterType","ST_FILTER_",NULL};
18: const char *STFilterDampings[] = {"NONE","JACKSON","LANCZOS","FEJER","STFilterDamping","ST_FILTER_DAMPING_",NULL};
20: static PetscErrorCode STFilterSetType_Private(ST st,STFilterType type)
21: {
22: ST_FILTER *ctx = (ST_FILTER*)st->data;
24: PetscFunctionBegin;
25: ctx->type = type;
26: switch(type) {
27: case ST_FILTER_FILTLAN:
28: PetscCall(STCreate_Filter_FILTLAN(st));
29: break;
30: case ST_FILTER_CHEBYSHEV:
31: PetscCall(STCreate_Filter_Chebyshev(st));
32: break;
33: default:
34: SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Invalid type");
35: }
36: PetscFunctionReturn(PETSC_SUCCESS);
37: }
39: /*
40: Operator (filter):
41: Op P M
42: if nmat=1: p(A) NULL p(A)
43: */
44: static PetscErrorCode STComputeOperator_Filter(ST st)
45: {
46: ST_FILTER *ctx = (ST_FILTER*)st->data;
48: PetscFunctionBegin;
49: PetscCheck(st->nmat==1,PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Only implemented for standard eigenvalue problem");
50: PetscCheck(ctx->intb<PETSC_MAX_REAL || ctx->inta>PETSC_MIN_REAL,PetscObjectComm((PetscObject)st),PETSC_ERR_ORDER,"Must pass an interval with STFilterSetInterval()");
51: PetscCheck(ctx->right!=0.0 || ctx->left!=0.0,PetscObjectComm((PetscObject)st),PETSC_ERR_ORDER,"Must pass an approximate numerical range with STFilterSetRange()");
52: PetscCheck(ctx->left<=ctx->inta && ctx->right>=ctx->intb,PetscObjectComm((PetscObject)st),PETSC_ERR_USER_INPUT,"The requested interval [%g,%g] must be contained in the numerical range [%g,%g]",(double)ctx->inta,(double)ctx->intb,(double)ctx->left,(double)ctx->right);
54: if (!ctx->type) PetscCall(STFilterSetType_Private(st,ST_FILTER_FILTLAN)); /* default type */
55: if (!ctx->polyDegree) ctx->polyDegree = 100;
56: PetscCall(ctx->computeoperator(st,&st->T[0]));
57: st->M = st->T[0];
58: PetscCall(MatDestroy(&st->P));
59: PetscFunctionReturn(PETSC_SUCCESS);
60: }
62: static PetscErrorCode STSetFromOptions_Filter(ST st,PetscOptionItems PetscOptionsObject)
63: {
64: ST_FILTER *ctx = (ST_FILTER*)st->data;
65: PetscReal array[2]={0,0};
66: PetscInt k;
67: PetscBool flg;
68: STFilterType type;
69: STFilterDamping damping;
71: PetscFunctionBegin;
72: PetscOptionsHeadBegin(PetscOptionsObject,"ST Filter Options");
74: PetscCall(PetscOptionsEnum("-st_filter_type","How to construct the filter","STFilterSetType",STFilterTypes,(PetscEnum)ctx->type,(PetscEnum*)&type,&flg));
75: if (flg) PetscCall(STFilterSetType(st,type));
77: k = 2;
78: PetscCall(PetscOptionsRealArray("-st_filter_interval","Interval containing the desired eigenvalues (two real values separated with a comma without spaces)","STFilterSetInterval",array,&k,&flg));
79: if (flg) {
80: PetscCheck(k>1,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_SIZ,"Must pass two values in -st_filter_interval (comma-separated without spaces)");
81: PetscCall(STFilterSetInterval(st,array[0],array[1]));
82: }
83: k = 2;
84: PetscCall(PetscOptionsRealArray("-st_filter_range","Interval containing all eigenvalues (two real values separated with a comma without spaces)","STFilterSetRange",array,&k,&flg));
85: if (flg) {
86: PetscCheck(k>1,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_SIZ,"Must pass two values in -st_filter_range (comma-separated without spaces)");
87: PetscCall(STFilterSetRange(st,array[0],array[1]));
88: }
89: PetscCall(PetscOptionsInt("-st_filter_degree","Degree of filter polynomial","STFilterSetDegree",100,&k,&flg));
90: if (flg) PetscCall(STFilterSetDegree(st,k));
92: PetscCall(PetscOptionsEnum("-st_filter_damping","Type of damping","STFilterSetDamping",STFilterDampings,(PetscEnum)ctx->damping,(PetscEnum*)&damping,&flg));
93: if (flg) PetscCall(STFilterSetDamping(st,damping));
95: PetscOptionsHeadEnd();
96: PetscFunctionReturn(PETSC_SUCCESS);
97: }
99: static PetscErrorCode STReset_Filter(ST st)
100: {
101: ST_FILTER *ctx = (ST_FILTER*)st->data;
103: PetscFunctionBegin;
104: ctx->left = 0.0;
105: ctx->right = 0.0;
106: PetscCall(MatDestroy(&ctx->T));
107: PetscCall(MatDestroyMatrices(ctx->nW,&ctx->W));
108: ctx->nW = 0;
109: PetscFunctionReturn(PETSC_SUCCESS);
110: }
112: static PetscErrorCode STFilterSetType_Filter(ST st,STFilterType type)
113: {
114: ST_FILTER *ctx = (ST_FILTER*)st->data;
116: PetscFunctionBegin;
117: if (ctx->type != type) {
118: PetscCall(STReset_Filter(st));
119: PetscCall(STFilterSetType_Private(st,type));
120: st->state = ST_STATE_INITIAL;
121: st->opready = PETSC_FALSE;
122: ctx->filtch = PETSC_TRUE;
123: }
124: PetscFunctionReturn(PETSC_SUCCESS);
125: }
127: /*@
128: STFilterSetType - Sets the method to be used to build the polynomial filter.
130: Logically Collective
132: Input Parameters:
133: + st - the spectral transformation context
134: - type - the type of filter
136: Options Database Key:
137: . -st_filter_type <type> - set the type of filter
139: Level: intermediate
141: .seealso: [](ch:st), `STFILTER`, `STFilterGetType()`
142: @*/
143: PetscErrorCode STFilterSetType(ST st,STFilterType type)
144: {
145: PetscFunctionBegin;
148: PetscTryMethod(st,"STFilterSetType_C",(ST,STFilterType),(st,type));
149: PetscFunctionReturn(PETSC_SUCCESS);
150: }
152: static PetscErrorCode STFilterGetType_Filter(ST st,STFilterType *type)
153: {
154: ST_FILTER *ctx = (ST_FILTER*)st->data;
156: PetscFunctionBegin;
157: *type = ctx->type ;
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
161: /*@
162: STFilterGetType - Gets the method to be used to build the polynomial filter.
164: Not Collective
166: Input Parameter:
167: . st - the spectral transformation context
169: Output Parameter:
170: . type - the type of filter
172: Level: intermediate
174: .seealso: [](ch:st), `STFilterSetType()`
175: @*/
176: PetscErrorCode STFilterGetType(ST st,STFilterType *type)
177: {
178: PetscFunctionBegin;
180: PetscAssertPointer(type,2);
181: PetscUseMethod(st,"STFilterGetType_C",(ST,STFilterType*),(st,type));
182: PetscFunctionReturn(PETSC_SUCCESS);
183: }
185: static PetscErrorCode STFilterSetInterval_Filter(ST st,PetscReal inta,PetscReal intb)
186: {
187: ST_FILTER *ctx = (ST_FILTER*)st->data;
189: PetscFunctionBegin;
190: PetscCheck(inta<intb,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
191: if (ctx->inta != inta || ctx->intb != intb) {
192: ctx->inta = inta;
193: ctx->intb = intb;
194: st->state = ST_STATE_INITIAL;
195: st->opready = PETSC_FALSE;
196: ctx->filtch = PETSC_TRUE;
197: }
198: PetscFunctionReturn(PETSC_SUCCESS);
199: }
201: /*@
202: STFilterSetInterval - Defines the interval containing the desired eigenvalues.
204: Logically Collective
206: Input Parameters:
207: + st - the spectral transformation context
208: . inta - left end of the interval
209: - intb - right end of the interval
211: Options Database Key:
212: . -st_filter_interval <a,b> - set [a,b] as the interval of interest
214: Notes:
215: The filter will be configured to emphasize eigenvalues contained in the given
216: interval, and damp out eigenvalues outside it. If the interval is open, then
217: the filter is low- or high-pass, otherwise it is mid-pass.
219: Common usage is to set the interval in `EPS` with `EPSSetInterval()`.
221: The interval must be contained within the numerical range of the matrix, see
222: `STFilterSetRange()`.
224: Level: intermediate
226: .seealso: [](ch:st), `STFILTER`, `STFilterGetInterval()`, `STFilterSetRange()`, `EPSSetInterval()`
227: @*/
228: PetscErrorCode STFilterSetInterval(ST st,PetscReal inta,PetscReal intb)
229: {
230: PetscFunctionBegin;
234: PetscTryMethod(st,"STFilterSetInterval_C",(ST,PetscReal,PetscReal),(st,inta,intb));
235: PetscFunctionReturn(PETSC_SUCCESS);
236: }
238: static PetscErrorCode STFilterGetInterval_Filter(ST st,PetscReal *inta,PetscReal *intb)
239: {
240: ST_FILTER *ctx = (ST_FILTER*)st->data;
242: PetscFunctionBegin;
243: if (inta) *inta = ctx->inta;
244: if (intb) *intb = ctx->intb;
245: PetscFunctionReturn(PETSC_SUCCESS);
246: }
248: /*@
249: STFilterGetInterval - Gets the interval containing the desired eigenvalues.
251: Not Collective
253: Input Parameter:
254: . st - the spectral transformation context
256: Output Parameters:
257: + inta - left end of the interval
258: - intb - right end of the interval
260: Level: intermediate
262: .seealso: [](ch:st), `STFilterSetInterval()`
263: @*/
264: PetscErrorCode STFilterGetInterval(ST st,PetscReal *inta,PetscReal *intb)
265: {
266: PetscFunctionBegin;
268: PetscUseMethod(st,"STFilterGetInterval_C",(ST,PetscReal*,PetscReal*),(st,inta,intb));
269: PetscFunctionReturn(PETSC_SUCCESS);
270: }
272: static PetscErrorCode STFilterSetRange_Filter(ST st,PetscReal left,PetscReal right)
273: {
274: ST_FILTER *ctx = (ST_FILTER*)st->data;
276: PetscFunctionBegin;
277: PetscCheck(left<right,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be left<right");
278: if (ctx->left != left || ctx->right != right) {
279: ctx->left = left;
280: ctx->right = right;
281: st->state = ST_STATE_INITIAL;
282: st->opready = PETSC_FALSE;
283: ctx->filtch = PETSC_TRUE;
284: }
285: PetscFunctionReturn(PETSC_SUCCESS);
286: }
288: /*@
289: STFilterSetRange - Defines the numerical range (or field of values) of the matrix, that is,
290: the interval containing all eigenvalues.
292: Logically Collective
294: Input Parameters:
295: + st - the spectral transformation context
296: . left - left end of the spectral range
297: - right - right end of the spectral range
299: Options Database Key:
300: . -st_filter_range <a,b> - set [a,b] as the numerical range
302: Notes:
303: The filter will be most effective if the numerical range is tight, that is,
304: `left` and `right` are good approximations to the leftmost and rightmost
305: eigenvalues, respectively.
307: Level: intermediate
309: .seealso: [](ch:st), `STFILTER`, `STFilterGetRange()`, `STFilterSetInterval()`
310: @*/
311: PetscErrorCode STFilterSetRange(ST st,PetscReal left,PetscReal right)
312: {
313: PetscFunctionBegin;
317: PetscTryMethod(st,"STFilterSetRange_C",(ST,PetscReal,PetscReal),(st,left,right));
318: PetscFunctionReturn(PETSC_SUCCESS);
319: }
321: static PetscErrorCode STFilterGetRange_Filter(ST st,PetscReal *left,PetscReal *right)
322: {
323: ST_FILTER *ctx = (ST_FILTER*)st->data;
325: PetscFunctionBegin;
326: if (left) *left = ctx->left;
327: if (right) *right = ctx->right;
328: PetscFunctionReturn(PETSC_SUCCESS);
329: }
331: /*@
332: STFilterGetRange - Gets the interval containing all eigenvalues.
334: Not Collective
336: Input Parameter:
337: . st - the spectral transformation context
339: Output Parameters:
340: + left - left end of the spectral range
341: - right - right end of the spectral range
343: Level: intermediate
345: .seealso: [](ch:st), `STFilterSetRange()`
346: @*/
347: PetscErrorCode STFilterGetRange(ST st,PetscReal *left,PetscReal *right)
348: {
349: PetscFunctionBegin;
351: PetscUseMethod(st,"STFilterGetRange_C",(ST,PetscReal*,PetscReal*),(st,left,right));
352: PetscFunctionReturn(PETSC_SUCCESS);
353: }
355: static PetscErrorCode STFilterSetDegree_Filter(ST st,PetscInt deg)
356: {
357: ST_FILTER *ctx = (ST_FILTER*)st->data;
359: PetscFunctionBegin;
360: if (deg == PETSC_DEFAULT || deg == PETSC_DECIDE) {
361: ctx->polyDegree = 0;
362: st->state = ST_STATE_INITIAL;
363: st->opready = PETSC_FALSE;
364: ctx->filtch = PETSC_TRUE;
365: } else {
366: PetscCheck(deg>0,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of degree. Must be > 0");
367: if (ctx->polyDegree != deg) {
368: ctx->polyDegree = deg;
369: st->state = ST_STATE_INITIAL;
370: st->opready = PETSC_FALSE;
371: ctx->filtch = PETSC_TRUE;
372: }
373: }
374: PetscFunctionReturn(PETSC_SUCCESS);
375: }
377: /*@
378: STFilterSetDegree - Sets the degree of the filter polynomial.
380: Logically Collective
382: Input Parameters:
383: + st - the spectral transformation context
384: - deg - polynomial degree
386: Options Database Key:
387: . -st_filter_degree <deg> - sets the degree of the filter polynomial
389: Level: intermediate
391: .seealso: [](ch:st), `STFILTER`, `STFilterGetDegree()`
392: @*/
393: PetscErrorCode STFilterSetDegree(ST st,PetscInt deg)
394: {
395: PetscFunctionBegin;
398: PetscTryMethod(st,"STFilterSetDegree_C",(ST,PetscInt),(st,deg));
399: PetscFunctionReturn(PETSC_SUCCESS);
400: }
402: static PetscErrorCode STFilterGetDegree_Filter(ST st,PetscInt *deg)
403: {
404: ST_FILTER *ctx = (ST_FILTER*)st->data;
406: PetscFunctionBegin;
407: *deg = ctx->polyDegree;
408: PetscFunctionReturn(PETSC_SUCCESS);
409: }
411: /*@
412: STFilterGetDegree - Gets the degree of the filter polynomial.
414: Not Collective
416: Input Parameter:
417: . st - the spectral transformation context
419: Output Parameter:
420: . deg - polynomial degree
422: Level: intermediate
424: .seealso: [](ch:st), `STFilterSetDegree()`
425: @*/
426: PetscErrorCode STFilterGetDegree(ST st,PetscInt *deg)
427: {
428: PetscFunctionBegin;
430: PetscAssertPointer(deg,2);
431: PetscUseMethod(st,"STFilterGetDegree_C",(ST,PetscInt*),(st,deg));
432: PetscFunctionReturn(PETSC_SUCCESS);
433: }
435: static PetscErrorCode STFilterGetThreshold_Filter(ST st,PetscReal *gamma)
436: {
437: ST_FILTER *ctx = (ST_FILTER*)st->data;
439: PetscFunctionBegin;
440: if (ctx->getthreshold) PetscCall(ctx->getthreshold(st,gamma));
441: else *gamma = 0.5;
442: PetscFunctionReturn(PETSC_SUCCESS);
443: }
445: /*@
446: STFilterGetThreshold - Gets the threshold value $\gamma$ used to decide
447: whether eigenvalue approximations are inside or outside the wanted interval.
449: Not Collective
451: Input Parameter:
452: . st - the spectral transformation context
454: Output Parameter:
455: . gamma - the threshold value
457: Note:
458: An eigenvalue $\lambda$ is considered to be inside the wanted interval
459: if $|p(\lambda)|\geq\gamma$, and outside otherwise.
461: Level: developer
463: .seealso: [](ch:st), `STFILTER`, `STFilterGetRange()`
464: @*/
465: PetscErrorCode STFilterGetThreshold(ST st,PetscReal *gamma)
466: {
467: PetscFunctionBegin;
469: PetscAssertPointer(gamma,2);
470: PetscUseMethod(st,"STFilterGetThreshold_C",(ST,PetscReal*),(st,gamma));
471: PetscFunctionReturn(PETSC_SUCCESS);
472: }
474: static PetscErrorCode STFilterSetDamping_Filter(ST st,STFilterDamping damping)
475: {
476: ST_FILTER *ctx = (ST_FILTER*)st->data;
478: PetscFunctionBegin;
479: if (ctx->damping != damping) {
480: ctx->damping = damping;
481: st->state = ST_STATE_INITIAL;
482: st->opready = PETSC_FALSE;
483: ctx->filtch = PETSC_TRUE;
484: }
485: PetscFunctionReturn(PETSC_SUCCESS);
486: }
488: /*@
489: STFilterSetDamping - Sets the type of damping to be used in the polynomial filter.
491: Logically Collective
493: Input Parameters:
494: + st - the spectral transformation context
495: - damping - the type of damping
497: Options Database Key:
498: . -st_filter_damping <damping> - sets the type of damping
500: Note:
501: Only used in `ST_FILTER_CHEBYSHEV` filters.
503: Level: advanced
505: .seealso: [](ch:st), `STFILTER`, `STFilterGetDamping()`
506: @*/
507: PetscErrorCode STFilterSetDamping(ST st,STFilterDamping damping)
508: {
509: PetscFunctionBegin;
512: PetscTryMethod(st,"STFilterSetDamping_C",(ST,STFilterDamping),(st,damping));
513: PetscFunctionReturn(PETSC_SUCCESS);
514: }
516: static PetscErrorCode STFilterGetDamping_Filter(ST st,STFilterDamping *damping)
517: {
518: ST_FILTER *ctx = (ST_FILTER*)st->data;
520: PetscFunctionBegin;
521: *damping = ctx->damping;
522: PetscFunctionReturn(PETSC_SUCCESS);
523: }
525: /*@
526: STFilterGetDamping - Gets the type of damping used in the polynomial filter.
528: Not Collective
530: Input Parameter:
531: . st - the spectral transformation context
533: Output Parameter:
534: . damping - the type of damping
536: Level: advanced
538: .seealso: [](ch:st), `STFilterSetDamping()`
539: @*/
540: PetscErrorCode STFilterGetDamping(ST st,STFilterDamping *damping)
541: {
542: PetscFunctionBegin;
544: PetscAssertPointer(damping,2);
545: PetscUseMethod(st,"STFilterGetDamping_C",(ST,STFilterDamping*),(st,damping));
546: PetscFunctionReturn(PETSC_SUCCESS);
547: }
549: static PetscErrorCode STView_Filter(ST st,PetscViewer viewer)
550: {
551: ST_FILTER *ctx = (ST_FILTER*)st->data;
552: PetscReal gamma;
553: PetscBool isascii;
555: PetscFunctionBegin;
556: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
557: if (isascii) {
558: PetscCall(PetscViewerASCIIPrintf(viewer," filter type: %s\n",STFilterTypes[ctx->type]));
559: PetscCall(PetscViewerASCIIPrintf(viewer," interval of desired eigenvalues: [%g,%g]\n",(double)ctx->inta,(double)ctx->intb));
560: PetscCall(PetscViewerASCIIPrintf(viewer," numerical range: [%g,%g]\n",(double)ctx->left,(double)ctx->right));
561: PetscCall(PetscViewerASCIIPrintf(viewer," degree of filter polynomial: %" PetscInt_FMT "\n",ctx->polyDegree));
562: if (ctx->damping && ctx->type==ST_FILTER_CHEBYSHEV) PetscCall(PetscViewerASCIIPrintf(viewer," type of damping = %s\n",STFilterDampings[ctx->damping]));
563: if (st->state>=ST_STATE_SETUP && ctx->getthreshold) {
564: PetscCall(STFilterGetThreshold(st,&gamma));
565: PetscCall(PetscViewerASCIIPrintf(viewer," limit to accept eigenvalues: gamma=%g\n",(double)gamma));
566: }
567: }
568: PetscFunctionReturn(PETSC_SUCCESS);
569: }
571: static PetscErrorCode STDestroy_Filter(ST st)
572: {
573: ST_FILTER *ctx = (ST_FILTER*)st->data;
575: PetscFunctionBegin;
576: if (ctx->destroy) PetscCall(ctx->destroy(st));
577: PetscCall(PetscFree(st->data));
578: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetType_C",NULL));
579: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetType_C",NULL));
580: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetInterval_C",NULL));
581: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetInterval_C",NULL));
582: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetRange_C",NULL));
583: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetRange_C",NULL));
584: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDegree_C",NULL));
585: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDegree_C",NULL));
586: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetThreshold_C",NULL));
587: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDamping_C",NULL));
588: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDamping_C",NULL));
589: PetscFunctionReturn(PETSC_SUCCESS);
590: }
592: /*MC
593: STFILTER - STFILTER = "filter" - A special type of `ST` that represents
594: a polynomial filter.
596: Level: beginner
598: Notes:
599: In standard eigenvalue problems, when the eigenvalues of interest are interior
600: to the spectrum and we want to avoid the high cost associated with the matrix
601: factorization of the shift-and-invert spectral transformation, an alternative
602: is to build a high-order polynomial such that $p(A)$ enhances the wanted
603: eigenvalues and filters out the unwanted ones.
605: The definition of the filter is done with functions such as `STFilterSetDegree()`,
606: `STFilterSetInterval()` or `STFilterSetType()`.
608: .seealso: [](ch:st), `ST`, `STType`, `STSetType()`, `STSetMatrices()`, `STFilterSetDegree()`, `STFilterSetInterval()`, `STFilterSetType()`
609: M*/
611: SLEPC_EXTERN PetscErrorCode STCreate_Filter(ST st)
612: {
613: ST_FILTER *ctx;
615: PetscFunctionBegin;
616: PetscCall(PetscNew(&ctx));
617: st->data = (void*)ctx;
619: st->usesksp = PETSC_FALSE;
621: ctx->type = (STFilterType)0;
622: ctx->inta = PETSC_MIN_REAL;
623: ctx->intb = PETSC_MAX_REAL;
624: ctx->left = 0.0;
625: ctx->right = 0.0;
626: ctx->polyDegree = 0;
628: st->ops->apply = STApply_Generic;
629: st->ops->computeoperator = STComputeOperator_Filter;
630: st->ops->setfromoptions = STSetFromOptions_Filter;
631: st->ops->destroy = STDestroy_Filter;
632: st->ops->reset = STReset_Filter;
633: st->ops->view = STView_Filter;
635: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetType_C",STFilterSetType_Filter));
636: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetType_C",STFilterGetType_Filter));
637: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetInterval_C",STFilterSetInterval_Filter));
638: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetInterval_C",STFilterGetInterval_Filter));
639: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetRange_C",STFilterSetRange_Filter));
640: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetRange_C",STFilterGetRange_Filter));
641: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDegree_C",STFilterSetDegree_Filter));
642: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDegree_C",STFilterGetDegree_Filter));
643: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetThreshold_C",STFilterGetThreshold_Filter));
644: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDamping_C",STFilterSetDamping_Filter));
645: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDamping_C",STFilterGetDamping_Filter));
646: PetscFunctionReturn(PETSC_SUCCESS);
647: }