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: `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: `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: `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: `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 interval
297: -  right - right end of the interval

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, left and right
304:    are good approximations to the leftmost and rightmost eigenvalues, respectively.

306:    Level: intermediate

308: .seealso: `STFilterGetRange()`, `STFilterSetInterval()`
309: @*/
310: PetscErrorCode STFilterSetRange(ST st,PetscReal left,PetscReal right)
311: {
312:   PetscFunctionBegin;
316:   PetscTryMethod(st,"STFilterSetRange_C",(ST,PetscReal,PetscReal),(st,left,right));
317:   PetscFunctionReturn(PETSC_SUCCESS);
318: }

320: static PetscErrorCode STFilterGetRange_Filter(ST st,PetscReal *left,PetscReal *right)
321: {
322:   ST_FILTER *ctx = (ST_FILTER*)st->data;

324:   PetscFunctionBegin;
325:   if (left)  *left  = ctx->left;
326:   if (right) *right = ctx->right;
327:   PetscFunctionReturn(PETSC_SUCCESS);
328: }

330: /*@
331:    STFilterGetRange - Gets the interval containing all eigenvalues.

333:    Not Collective

335:    Input Parameter:
336: .  st  - the spectral transformation context

338:    Output Parameters:
339: +  left  - left end of the interval
340: -  right - right end of the interval

342:    Level: intermediate

344: .seealso: `STFilterSetRange()`
345: @*/
346: PetscErrorCode STFilterGetRange(ST st,PetscReal *left,PetscReal *right)
347: {
348:   PetscFunctionBegin;
350:   PetscUseMethod(st,"STFilterGetRange_C",(ST,PetscReal*,PetscReal*),(st,left,right));
351:   PetscFunctionReturn(PETSC_SUCCESS);
352: }

354: static PetscErrorCode STFilterSetDegree_Filter(ST st,PetscInt deg)
355: {
356:   ST_FILTER *ctx = (ST_FILTER*)st->data;

358:   PetscFunctionBegin;
359:   if (deg == PETSC_DEFAULT || deg == PETSC_DECIDE) {
360:     ctx->polyDegree = 0;
361:     st->state       = ST_STATE_INITIAL;
362:     st->opready     = PETSC_FALSE;
363:     ctx->filtch     = PETSC_TRUE;
364:   } else {
365:     PetscCheck(deg>0,PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of degree. Must be > 0");
366:     if (ctx->polyDegree != deg) {
367:       ctx->polyDegree = deg;
368:       st->state       = ST_STATE_INITIAL;
369:       st->opready     = PETSC_FALSE;
370:       ctx->filtch     = PETSC_TRUE;
371:     }
372:   }
373:   PetscFunctionReturn(PETSC_SUCCESS);
374: }

376: /*@
377:    STFilterSetDegree - Sets the degree of the filter polynomial.

379:    Logically Collective

381:    Input Parameters:
382: +  st  - the spectral transformation context
383: -  deg - polynomial degree

385:    Options Database Key:
386: .  -st_filter_degree <deg> - sets the degree of the filter polynomial

388:    Level: intermediate

390: .seealso: `STFilterGetDegree()`
391: @*/
392: PetscErrorCode STFilterSetDegree(ST st,PetscInt deg)
393: {
394:   PetscFunctionBegin;
397:   PetscTryMethod(st,"STFilterSetDegree_C",(ST,PetscInt),(st,deg));
398:   PetscFunctionReturn(PETSC_SUCCESS);
399: }

401: static PetscErrorCode STFilterGetDegree_Filter(ST st,PetscInt *deg)
402: {
403:   ST_FILTER *ctx = (ST_FILTER*)st->data;

405:   PetscFunctionBegin;
406:   *deg = ctx->polyDegree;
407:   PetscFunctionReturn(PETSC_SUCCESS);
408: }

410: /*@
411:    STFilterGetDegree - Gets the degree of the filter polynomial.

413:    Not Collective

415:    Input Parameter:
416: .  st  - the spectral transformation context

418:    Output Parameter:
419: .  deg - polynomial degree

421:    Level: intermediate

423: .seealso: `STFilterSetDegree()`
424: @*/
425: PetscErrorCode STFilterGetDegree(ST st,PetscInt *deg)
426: {
427:   PetscFunctionBegin;
429:   PetscAssertPointer(deg,2);
430:   PetscUseMethod(st,"STFilterGetDegree_C",(ST,PetscInt*),(st,deg));
431:   PetscFunctionReturn(PETSC_SUCCESS);
432: }

434: static PetscErrorCode STFilterGetThreshold_Filter(ST st,PetscReal *gamma)
435: {
436:   ST_FILTER *ctx = (ST_FILTER*)st->data;

438:   PetscFunctionBegin;
439:   if (ctx->getthreshold) PetscCall(ctx->getthreshold(st,gamma));
440:   else *gamma = 0.5;
441:   PetscFunctionReturn(PETSC_SUCCESS);
442: }

444: /*@
445:    STFilterGetThreshold - Gets the threshold value gamma such that rho(lambda)>=gamma for
446:    eigenvalues lambda inside the wanted interval and rho(lambda)<gamma for those outside.

448:    Not Collective

450:    Input Parameter:
451: .  st  - the spectral transformation context

453:    Output Parameter:
454: .  gamma - the threshold value

456:    Level: developer

458: .seealso: `STFilterGetRange()`
459: @*/
460: PetscErrorCode STFilterGetThreshold(ST st,PetscReal *gamma)
461: {
462:   PetscFunctionBegin;
464:   PetscAssertPointer(gamma,2);
465:   PetscUseMethod(st,"STFilterGetThreshold_C",(ST,PetscReal*),(st,gamma));
466:   PetscFunctionReturn(PETSC_SUCCESS);
467: }

469: static PetscErrorCode STFilterSetDamping_Filter(ST st,STFilterDamping damping)
470: {
471:   ST_FILTER *ctx = (ST_FILTER*)st->data;

473:   PetscFunctionBegin;
474:   if (ctx->damping != damping) {
475:     ctx->damping = damping;
476:     st->state    = ST_STATE_INITIAL;
477:     st->opready  = PETSC_FALSE;
478:     ctx->filtch  = PETSC_TRUE;
479:   }
480:   PetscFunctionReturn(PETSC_SUCCESS);
481: }

483: /*@
484:    STFilterSetDamping - Sets the type of damping to be used in the polynomial filter.

486:    Logically Collective

488:    Input Parameters:
489: +  st      - the spectral transformation context
490: -  damping - the type of damping

492:    Options Database Key:
493: .  -st_filter_damping <damping> - sets the type of damping

495:    Note:
496:    Only used in Chebyshev filters.

498:    Level: advanced

500: .seealso: `STFilterGetDamping()`
501: @*/
502: PetscErrorCode STFilterSetDamping(ST st,STFilterDamping damping)
503: {
504:   PetscFunctionBegin;
507:   PetscTryMethod(st,"STFilterSetDamping_C",(ST,STFilterDamping),(st,damping));
508:   PetscFunctionReturn(PETSC_SUCCESS);
509: }

511: static PetscErrorCode STFilterGetDamping_Filter(ST st,STFilterDamping *damping)
512: {
513:   ST_FILTER *ctx = (ST_FILTER*)st->data;

515:   PetscFunctionBegin;
516:   *damping = ctx->damping;
517:   PetscFunctionReturn(PETSC_SUCCESS);
518: }

520: /*@
521:    STFilterGetDamping - Gets the type of damping used in the polynomial filter.

523:    Not Collective

525:    Input Parameter:
526: .  st  - the spectral transformation context

528:    Output Parameter:
529: .  damping - the type of damping

531:    Level: advanced

533: .seealso: `STFilterSetDamping()`
534: @*/
535: PetscErrorCode STFilterGetDamping(ST st,STFilterDamping *damping)
536: {
537:   PetscFunctionBegin;
539:   PetscAssertPointer(damping,2);
540:   PetscUseMethod(st,"STFilterGetDamping_C",(ST,STFilterDamping*),(st,damping));
541:   PetscFunctionReturn(PETSC_SUCCESS);
542: }

544: static PetscErrorCode STView_Filter(ST st,PetscViewer viewer)
545: {
546:   ST_FILTER *ctx = (ST_FILTER*)st->data;
547:   PetscReal gamma;
548:   PetscBool isascii;

550:   PetscFunctionBegin;
551:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
552:   if (isascii) {
553:     PetscCall(PetscViewerASCIIPrintf(viewer,"  filter type: %s\n",STFilterTypes[ctx->type]));
554:     PetscCall(PetscViewerASCIIPrintf(viewer,"  interval of desired eigenvalues: [%g,%g]\n",(double)ctx->inta,(double)ctx->intb));
555:     PetscCall(PetscViewerASCIIPrintf(viewer,"  numerical range: [%g,%g]\n",(double)ctx->left,(double)ctx->right));
556:     PetscCall(PetscViewerASCIIPrintf(viewer,"  degree of filter polynomial: %" PetscInt_FMT "\n",ctx->polyDegree));
557:     if (ctx->damping && ctx->type==ST_FILTER_CHEBYSHEV) PetscCall(PetscViewerASCIIPrintf(viewer,"  type of damping = %s\n",STFilterDampings[ctx->damping]));
558:     if (st->state>=ST_STATE_SETUP && ctx->getthreshold) {
559:       PetscCall(STFilterGetThreshold(st,&gamma));
560:       PetscCall(PetscViewerASCIIPrintf(viewer,"  limit to accept eigenvalues: gamma=%g\n",(double)gamma));
561:     }
562:   }
563:   PetscFunctionReturn(PETSC_SUCCESS);
564: }

566: static PetscErrorCode STDestroy_Filter(ST st)
567: {
568:   ST_FILTER *ctx = (ST_FILTER*)st->data;

570:   PetscFunctionBegin;
571:   if (ctx->destroy) PetscCall(ctx->destroy(st));
572:   PetscCall(PetscFree(st->data));
573:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetType_C",NULL));
574:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetType_C",NULL));
575:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetInterval_C",NULL));
576:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetInterval_C",NULL));
577:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetRange_C",NULL));
578:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetRange_C",NULL));
579:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDegree_C",NULL));
580:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDegree_C",NULL));
581:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetThreshold_C",NULL));
582:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDamping_C",NULL));
583:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDamping_C",NULL));
584:   PetscFunctionReturn(PETSC_SUCCESS);
585: }

587: SLEPC_EXTERN PetscErrorCode STCreate_Filter(ST st)
588: {
589:   ST_FILTER   *ctx;

591:   PetscFunctionBegin;
592:   PetscCall(PetscNew(&ctx));
593:   st->data = (void*)ctx;

595:   st->usesksp = PETSC_FALSE;

597:   ctx->type               = (STFilterType)0;
598:   ctx->inta               = PETSC_MIN_REAL;
599:   ctx->intb               = PETSC_MAX_REAL;
600:   ctx->left               = 0.0;
601:   ctx->right              = 0.0;
602:   ctx->polyDegree         = 0;

604:   st->ops->apply           = STApply_Generic;
605:   st->ops->computeoperator = STComputeOperator_Filter;
606:   st->ops->setfromoptions  = STSetFromOptions_Filter;
607:   st->ops->destroy         = STDestroy_Filter;
608:   st->ops->reset           = STReset_Filter;
609:   st->ops->view            = STView_Filter;

611:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetType_C",STFilterSetType_Filter));
612:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetType_C",STFilterGetType_Filter));
613:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetInterval_C",STFilterSetInterval_Filter));
614:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetInterval_C",STFilterGetInterval_Filter));
615:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetRange_C",STFilterSetRange_Filter));
616:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetRange_C",STFilterGetRange_Filter));
617:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDegree_C",STFilterSetDegree_Filter));
618:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDegree_C",STFilterGetDegree_Filter));
619:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetThreshold_C",STFilterGetThreshold_Filter));
620:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterSetDamping_C",STFilterSetDamping_Filter));
621:   PetscCall(PetscObjectComposeFunction((PetscObject)st,"STFilterGetDamping_C",STFilterGetDamping_Filter));
622:   PetscFunctionReturn(PETSC_SUCCESS);
623: }