Actual source code: dsbasic.c

slepc-3.7.2 2016-07-19
Report Typos and Errors
  1: /*
  2:    Basic DS routines

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.

 10:    SLEPc is free software: you can redistribute it and/or modify it under  the
 11:    terms of version 3 of the GNU Lesser General Public License as published by
 12:    the Free Software Foundation.

 14:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 15:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 16:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 17:    more details.

 19:    You  should have received a copy of the GNU Lesser General  Public  License
 20:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 21:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 22: */

 24: #include <slepc/private/dsimpl.h>      /*I "slepcds.h" I*/

 26: PetscFunctionList DSList = 0;
 27: PetscBool         DSRegisterAllCalled = PETSC_FALSE;
 28: PetscClassId      DS_CLASSID = 0;
 29: PetscLogEvent     DS_Solve = 0,DS_Vectors = 0,DS_Other = 0;
 30: static PetscBool  DSPackageInitialized = PETSC_FALSE;
 31: const char        *DSMatName[DS_NUM_MAT] = {"A","B","C","T","D","Q","Z","X","Y","U","VT","W","E0","E1","E2","E3","E4","E5","E6","E7","E8","E9"};
 32: DSMatType         DSMatExtra[DS_NUM_EXTRA] = {DS_MAT_E0,DS_MAT_E1,DS_MAT_E2,DS_MAT_E3,DS_MAT_E4,DS_MAT_E5,DS_MAT_E6,DS_MAT_E7,DS_MAT_E8,DS_MAT_E9};

 36: /*@C
 37:    DSFinalizePackage - This function destroys everything in the SLEPc interface
 38:    to the DS package. It is called from SlepcFinalize().

 40:    Level: developer

 42: .seealso: SlepcFinalize()
 43: @*/
 44: PetscErrorCode DSFinalizePackage(void)
 45: {

 49:   PetscFunctionListDestroy(&DSList);
 50:   DSPackageInitialized = PETSC_FALSE;
 51:   DSRegisterAllCalled  = PETSC_FALSE;
 52:   return(0);
 53: }

 57: /*@C
 58:   DSInitializePackage - This function initializes everything in the DS package.
 59:   It is called from PetscDLLibraryRegister() when using dynamic libraries, and
 60:   on the first call to DSCreate() when using static libraries.

 62:   Level: developer

 64: .seealso: SlepcInitialize()
 65: @*/
 66: PetscErrorCode DSInitializePackage()
 67: {
 68:   char             logList[256];
 69:   char             *className;
 70:   PetscBool        opt;
 71:   PetscErrorCode   ierr;

 74:   if (DSPackageInitialized) return(0);
 75:   DSPackageInitialized = PETSC_TRUE;
 76:   /* Register Classes */
 77:   PetscClassIdRegister("Direct Solver",&DS_CLASSID);
 78:   /* Register Constructors */
 79:   DSRegisterAll();
 80:   /* Register Events */
 81:   PetscLogEventRegister("DSSolve",DS_CLASSID,&DS_Solve);
 82:   PetscLogEventRegister("DSVectors",DS_CLASSID,&DS_Vectors);
 83:   PetscLogEventRegister("DSOther",DS_CLASSID,&DS_Other);
 84:   /* Process info exclusions */
 85:   PetscOptionsGetString(NULL,NULL,"-info_exclude",logList,256,&opt);
 86:   if (opt) {
 87:     PetscStrstr(logList,"ds",&className);
 88:     if (className) {
 89:       PetscInfoDeactivateClass(DS_CLASSID);
 90:     }
 91:   }
 92:   /* Process summary exclusions */
 93:   PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,256,&opt);
 94:   if (opt) {
 95:     PetscStrstr(logList,"ds",&className);
 96:     if (className) {
 97:       PetscLogEventDeactivateClass(DS_CLASSID);
 98:     }
 99:   }
100:   PetscRegisterFinalize(DSFinalizePackage);
101:   return(0);
102: }

106: /*@
107:    DSCreate - Creates a DS context.

109:    Collective on MPI_Comm

111:    Input Parameter:
112: .  comm - MPI communicator

114:    Output Parameter:
115: .  newds - location to put the DS context

117:    Level: beginner

119:    Note:
120:    DS objects are not intended for normal users but only for
121:    advanced user that for instance implement their own solvers.

123: .seealso: DSDestroy(), DS
124: @*/
125: PetscErrorCode DSCreate(MPI_Comm comm,DS *newds)
126: {
127:   DS             ds;
128:   PetscInt       i;

133:   *newds = 0;
134:   DSInitializePackage();
135:   SlepcHeaderCreate(ds,DS_CLASSID,"DS","Direct Solver (or Dense System)","DS",comm,DSDestroy,DSView);

137:   ds->state         = DS_STATE_RAW;
138:   ds->method        = 0;
139:   ds->compact       = PETSC_FALSE;
140:   ds->refined       = PETSC_FALSE;
141:   ds->extrarow      = PETSC_FALSE;
142:   ds->ld            = 0;
143:   ds->l             = 0;
144:   ds->n             = 0;
145:   ds->m             = 0;
146:   ds->k             = 0;
147:   ds->t             = 0;
148:   ds->bs            = 1;
149:   ds->sc            = NULL;

151:   for (i=0;i<DS_NUM_MAT;i++) {
152:     ds->mat[i]      = NULL;
153:     ds->rmat[i]     = NULL;
154:     ds->omat[i]     = NULL;
155:   }
156:   ds->perm          = NULL;
157:   ds->data          = NULL;
158:   ds->work          = NULL;
159:   ds->rwork         = NULL;
160:   ds->iwork         = NULL;
161:   ds->lwork         = 0;
162:   ds->lrwork        = 0;
163:   ds->liwork        = 0;

165:   *newds = ds;
166:   return(0);
167: }

171: /*@C
172:    DSSetOptionsPrefix - Sets the prefix used for searching for all
173:    DS options in the database.

175:    Logically Collective on DS

177:    Input Parameters:
178: +  ds - the direct solver context
179: -  prefix - the prefix string to prepend to all DS option requests

181:    Notes:
182:    A hyphen (-) must NOT be given at the beginning of the prefix name.
183:    The first character of all runtime options is AUTOMATICALLY the
184:    hyphen.

186:    Level: advanced

188: .seealso: DSAppendOptionsPrefix()
189: @*/
190: PetscErrorCode DSSetOptionsPrefix(DS ds,const char *prefix)
191: {

196:   PetscObjectSetOptionsPrefix((PetscObject)ds,prefix);
197:   return(0);
198: }

202: /*@C
203:    DSAppendOptionsPrefix - Appends to the prefix used for searching for all
204:    DS options in the database.

206:    Logically Collective on DS

208:    Input Parameters:
209: +  ds - the direct solver context
210: -  prefix - the prefix string to prepend to all DS option requests

212:    Notes:
213:    A hyphen (-) must NOT be given at the beginning of the prefix name.
214:    The first character of all runtime options is AUTOMATICALLY the hyphen.

216:    Level: advanced

218: .seealso: DSSetOptionsPrefix()
219: @*/
220: PetscErrorCode DSAppendOptionsPrefix(DS ds,const char *prefix)
221: {

226:   PetscObjectAppendOptionsPrefix((PetscObject)ds,prefix);
227:   return(0);
228: }

232: /*@C
233:    DSGetOptionsPrefix - Gets the prefix used for searching for all
234:    DS options in the database.

236:    Not Collective

238:    Input Parameters:
239: .  ds - the direct solver context

241:    Output Parameters:
242: .  prefix - pointer to the prefix string used is returned

244:    Note:
245:    On the Fortran side, the user should pass in a string 'prefix' of
246:    sufficient length to hold the prefix.

248:    Level: advanced

250: .seealso: DSSetOptionsPrefix(), DSAppendOptionsPrefix()
251: @*/
252: PetscErrorCode DSGetOptionsPrefix(DS ds,const char *prefix[])
253: {

259:   PetscObjectGetOptionsPrefix((PetscObject)ds,prefix);
260:   return(0);
261: }

265: /*@C
266:    DSSetType - Selects the type for the DS object.

268:    Logically Collective on DS

270:    Input Parameter:
271: +  ds   - the direct solver context
272: -  type - a known type

274:    Level: intermediate

276: .seealso: DSGetType()
277: @*/
278: PetscErrorCode DSSetType(DS ds,DSType type)
279: {
280:   PetscErrorCode ierr,(*r)(DS);
281:   PetscBool      match;


287:   PetscObjectTypeCompare((PetscObject)ds,type,&match);
288:   if (match) return(0);

290:    PetscFunctionListFind(DSList,type,&r);
291:   if (!r) SETERRQ1(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested DS type %s",type);

293:   PetscMemzero(ds->ops,sizeof(struct _DSOps));

295:   PetscObjectChangeTypeName((PetscObject)ds,type);
296:   (*r)(ds);
297:   return(0);
298: }

302: /*@C
303:    DSGetType - Gets the DS type name (as a string) from the DS context.

305:    Not Collective

307:    Input Parameter:
308: .  ds - the direct solver context

310:    Output Parameter:
311: .  name - name of the direct solver

313:    Level: intermediate

315: .seealso: DSSetType()
316: @*/
317: PetscErrorCode DSGetType(DS ds,DSType *type)
318: {
322:   *type = ((PetscObject)ds)->type_name;
323:   return(0);
324: }

328: /*@
329:    DSSetMethod - Selects the method to be used to solve the problem.

331:    Logically Collective on DS

333:    Input Parameter:
334: +  ds   - the direct solver context
335: -  meth - an index indentifying the method

337:    Level: intermediate

339: .seealso: DSGetMethod()
340: @*/
341: PetscErrorCode DSSetMethod(DS ds,PetscInt meth)
342: {
346:   if (meth<0) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"The method must be a non-negative integer");
347:   if (meth>DS_MAX_SOLVE) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"Too large value for the method");
348:   ds->method = meth;
349:   return(0);
350: }

354: /*@
355:    DSGetMethod - Gets the method currently used in the DS.

357:    Not Collective

359:    Input Parameter:
360: .  ds - the direct solver context

362:    Output Parameter:
363: .  meth - identifier of the method

365:    Level: intermediate

367: .seealso: DSSetMethod()
368: @*/
369: PetscErrorCode DSGetMethod(DS ds,PetscInt *meth)
370: {
374:   *meth = ds->method;
375:   return(0);
376: }

380: /*@
381:    DSSetCompact - Switch to compact storage of matrices.

383:    Logically Collective on DS

385:    Input Parameter:
386: +  ds   - the direct solver context
387: -  comp - a boolean flag

389:    Notes:
390:    Compact storage is used in some DS types such as DSHEP when the matrix
391:    is tridiagonal. This flag can be used to indicate whether the user
392:    provides the matrix entries via the compact form (the tridiagonal DS_MAT_T)
393:    or the non-compact one (DS_MAT_A).

395:    The default is PETSC_FALSE.

397:    Level: advanced

399: .seealso: DSGetCompact()
400: @*/
401: PetscErrorCode DSSetCompact(DS ds,PetscBool comp)
402: {
406:   ds->compact = comp;
407:   return(0);
408: }

412: /*@
413:    DSGetCompact - Gets the compact storage flag.

415:    Not Collective

417:    Input Parameter:
418: .  ds - the direct solver context

420:    Output Parameter:
421: .  comp - the flag

423:    Level: advanced

425: .seealso: DSSetCompact()
426: @*/
427: PetscErrorCode DSGetCompact(DS ds,PetscBool *comp)
428: {
432:   *comp = ds->compact;
433:   return(0);
434: }

438: /*@
439:    DSSetExtraRow - Sets a flag to indicate that the matrix has one extra
440:    row.

442:    Logically Collective on DS

444:    Input Parameter:
445: +  ds  - the direct solver context
446: -  ext - a boolean flag

448:    Notes:
449:    In Krylov methods it is useful that the matrix representing the direct solver
450:    has one extra row, i.e., has dimension (n+1) x n. If this flag is activated, all
451:    transformations applied to the right of the matrix also affect this additional
452:    row. In that case, (n+1) must be less or equal than the leading dimension.

454:    The default is PETSC_FALSE.

456:    Level: advanced

458: .seealso: DSSolve(), DSAllocate(), DSGetExtraRow()
459: @*/
460: PetscErrorCode DSSetExtraRow(DS ds,PetscBool ext)
461: {
465:   if (ds->n>0 && ds->n==ds->ld) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ORDER,"Cannot set extra row after setting n=ld");
466:   ds->extrarow = ext;
467:   return(0);
468: }

472: /*@
473:    DSGetExtraRow - Gets the extra row flag.

475:    Not Collective

477:    Input Parameter:
478: .  ds - the direct solver context

480:    Output Parameter:
481: .  ext - the flag

483:    Level: advanced

485: .seealso: DSSetExtraRow()
486: @*/
487: PetscErrorCode DSGetExtraRow(DS ds,PetscBool *ext)
488: {
492:   *ext = ds->extrarow;
493:   return(0);
494: }

498: /*@
499:    DSSetRefined - Sets a flag to indicate that refined vectors must be
500:    computed.

502:    Logically Collective on DS

504:    Input Parameter:
505: +  ds  - the direct solver context
506: -  ref - a boolean flag

508:    Notes:
509:    Normally the vectors returned in DS_MAT_X are eigenvectors of the
510:    projected matrix. With this flag activated, DSVectors() will return
511:    the right singular vector of the smallest singular value of matrix
512:    \tilde{A}-theta*I, where \tilde{A} is the extended (n+1)xn matrix
513:    and theta is the Ritz value. This is used in the refined Ritz
514:    approximation.

516:    The default is PETSC_FALSE.

518:    Level: advanced

520: .seealso: DSVectors(), DSGetRefined()
521: @*/
522: PetscErrorCode DSSetRefined(DS ds,PetscBool ref)
523: {
527:   ds->refined = ref;
528:   return(0);
529: }

533: /*@
534:    DSGetRefined - Gets the refined vectors flag.

536:    Not Collective

538:    Input Parameter:
539: .  ds - the direct solver context

541:    Output Parameter:
542: .  ref - the flag

544:    Level: advanced

546: .seealso: DSSetRefined()
547: @*/
548: PetscErrorCode DSGetRefined(DS ds,PetscBool *ref)
549: {
553:   *ref = ds->refined;
554:   return(0);
555: }

559: /*@
560:    DSSetBlockSize - Sets the block size.

562:    Logically Collective on DS

564:    Input Parameter:
565: +  ds - the direct solver context
566: -  bs - the block size

568:    Level: intermediate

570: .seealso: DSGetBlockSize()
571: @*/
572: PetscErrorCode DSSetBlockSize(DS ds,PetscInt bs)
573: {
577:   if (bs<1) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"The block size must be at least one");
578:   ds->bs = bs;
579:   return(0);
580: }

584: /*@
585:    DSGetBlockSize - Gets the block size.

587:    Not Collective

589:    Input Parameter:
590: .  ds - the direct solver context

592:    Output Parameter:
593: .  bs - block size

595:    Level: intermediate

597: .seealso: DSSetBlockSize()
598: @*/
599: PetscErrorCode DSGetBlockSize(DS ds,PetscInt *bs)
600: {
604:   *bs = ds->bs;
605:   return(0);
606: }

610: /*@C
611:    DSSetSlepcSC - Sets the sorting criterion context.

613:    Not Collective

615:    Input Parameters:
616: +  ds - the direct solver context
617: -  sc - a pointer to the sorting criterion context

619:    Level: developer

621: .seealso: DSGetSlepcSC(), DSSort()
622: @*/
623: PetscErrorCode DSSetSlepcSC(DS ds,SlepcSC sc)
624: {

630:   if (ds->sc) {
631:     PetscFree(ds->sc);
632:   }
633:   ds->sc = sc;
634:   return(0);
635: }

639: /*@C
640:    DSGetSlepcSC - Gets the sorting criterion context.

642:    Not Collective

644:    Input Parameter:
645: .  ds - the direct solver context

647:    Output Parameters:
648: .  sc - a pointer to the sorting criterion context

650:    Level: developer

652: .seealso: DSSetSlepcSC(), DSSort()
653: @*/
654: PetscErrorCode DSGetSlepcSC(DS ds,SlepcSC *sc)
655: {

661:   if (!ds->sc) {
662:     PetscNewLog(ds,&ds->sc);
663:   }
664:   *sc = ds->sc;
665:   return(0);
666: }

670: /*@
671:    DSSetFromOptions - Sets DS options from the options database.

673:    Collective on DS

675:    Input Parameters:
676: .  ds - the direct solver context

678:    Notes:
679:    To see all options, run your program with the -help option.

681:    Level: beginner
682: @*/
683: PetscErrorCode DSSetFromOptions(DS ds)
684: {
686:   PetscInt       bs,meth;
687:   PetscBool      flag;

691:   DSRegisterAll();
692:   /* Set default type (we do not allow changing it with -ds_type) */
693:   if (!((PetscObject)ds)->type_name) {
694:     DSSetType(ds,DSNHEP);
695:   }
696:   PetscObjectOptionsBegin((PetscObject)ds);
697:     PetscOptionsInt("-ds_block_size","Block size for the dense system solver","DSSetBlockSize",ds->bs,&bs,&flag);
698:     if (flag) { DSSetBlockSize(ds,bs); }
699:     PetscOptionsInt("-ds_method","Method to be used for the dense system","DSSetMethod",ds->method,&meth,&flag);
700:     if (flag) { DSSetMethod(ds,meth); }
701:     PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)ds);
702:   PetscOptionsEnd();
703:   return(0);
704: }

708: /*@C
709:    DSView - Prints the DS data structure.

711:    Collective on DS

713:    Input Parameters:
714: +  ds - the direct solver context
715: -  viewer - optional visualization context

717:    Note:
718:    The available visualization contexts include
719: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
720: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
721:          output where only the first processor opens
722:          the file.  All other processors send their
723:          data to the first processor to print.

725:    The user can open an alternative visualization context with
726:    PetscViewerASCIIOpen() - output to a specified file.

728:    Level: beginner

730: .seealso: DSViewMat()
731: @*/
732: PetscErrorCode DSView(DS ds,PetscViewer viewer)
733: {
734:   PetscBool         isascii,issvd;
735:   const char        *state;
736:   PetscViewerFormat format;
737:   PetscErrorCode    ierr;

741:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ds));
744:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
745:   if (isascii) {
746:     PetscViewerGetFormat(viewer,&format);
747:     PetscObjectPrintClassNamePrefixType((PetscObject)ds,viewer);
748:     if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
749:       switch (ds->state) {
750:         case DS_STATE_RAW:          state = "raw"; break;
751:         case DS_STATE_INTERMEDIATE: state = "intermediate"; break;
752:         case DS_STATE_CONDENSED:    state = "condensed"; break;
753:         case DS_STATE_TRUNCATED:    state = "truncated"; break;
754:         default: SETERRQ(PetscObjectComm((PetscObject)ds),1,"Wrong value of ds->state");
755:       }
756:       PetscViewerASCIIPrintf(viewer,"  current state: %s\n",state);
757:       PetscObjectTypeCompare((PetscObject)ds,DSSVD,&issvd);
758:       if (issvd) {
759:         PetscViewerASCIIPrintf(viewer,"  dimensions: ld=%D, n=%D, m=%D, l=%D, k=%D",ds->ld,ds->n,ds->m,ds->l,ds->k);
760:       } else {
761:         PetscViewerASCIIPrintf(viewer,"  dimensions: ld=%D, n=%D, l=%D, k=%D",ds->ld,ds->n,ds->l,ds->k);
762:       }
763:       if (ds->state==DS_STATE_TRUNCATED) {
764:         PetscViewerASCIIPrintf(viewer,", t=%D\n",ds->t);
765:       } else {
766:         PetscViewerASCIIPrintf(viewer,"\n");
767:       }
768:       PetscViewerASCIIPrintf(viewer,"  flags:%s%s%s\n",ds->compact?" compact":"",ds->extrarow?" extrarow":"",ds->refined?" refined":"");
769:     }
770:     if (ds->ops->view) {
771:       PetscViewerASCIIPushTab(viewer);
772:       (*ds->ops->view)(ds,viewer);
773:       PetscViewerASCIIPopTab(viewer);
774:     }
775:   }
776:   return(0);
777: }

781: /*@
782:    DSAllocate - Allocates memory for internal storage or matrices in DS.

784:    Logically Collective on DS

786:    Input Parameters:
787: +  ds - the direct solver context
788: -  ld - leading dimension (maximum allowed dimension for the matrices, including
789:         the extra row if present)

791:    Level: intermediate

793: .seealso: DSGetLeadingDimension(), DSSetDimensions(), DSSetExtraRow()
794: @*/
795: PetscErrorCode DSAllocate(DS ds,PetscInt ld)
796: {

802:   if (ld<1) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"Leading dimension should be at least one");
803:   ds->ld = ld;
804:   (*ds->ops->allocate)(ds,ld);
805:   return(0);
806: }

810: /*@
811:    DSReset - Resets the DS context to the initial state.

813:    Collective on DS

815:    Input Parameter:
816: .  ds - the direct solver context

818:    Level: advanced

820: .seealso: DSDestroy()
821: @*/
822: PetscErrorCode DSReset(DS ds)
823: {
824:   PetscInt       i;

829:   ds->state    = DS_STATE_RAW;
830:   ds->compact  = PETSC_FALSE;
831:   ds->refined  = PETSC_FALSE;
832:   ds->extrarow = PETSC_FALSE;
833:   ds->ld       = 0;
834:   ds->l        = 0;
835:   ds->n        = 0;
836:   ds->m        = 0;
837:   ds->k        = 0;
838:   for (i=0;i<DS_NUM_MAT;i++) {
839:     PetscFree(ds->mat[i]);
840:     PetscFree(ds->rmat[i]);
841:     MatDestroy(&ds->omat[i]);
842:   }
843:   PetscFree(ds->perm);
844:   PetscFree(ds->work);
845:   PetscFree(ds->rwork);
846:   PetscFree(ds->iwork);
847:   ds->lwork         = 0;
848:   ds->lrwork        = 0;
849:   ds->liwork        = 0;
850:   return(0);
851: }

855: /*@
856:    DSDestroy - Destroys DS context that was created with DSCreate().

858:    Collective on DS

860:    Input Parameter:
861: .  ds - the direct solver context

863:    Level: beginner

865: .seealso: DSCreate()
866: @*/
867: PetscErrorCode DSDestroy(DS *ds)
868: {

872:   if (!*ds) return(0);
874:   if (--((PetscObject)(*ds))->refct > 0) { *ds = 0; return(0); }
875:   DSReset(*ds);
876:   if ((*ds)->ops->destroy) { (*(*ds)->ops->destroy)(*ds); }
877:   PetscFree((*ds)->sc);
878:   PetscHeaderDestroy(ds);
879:   return(0);
880: }

884: /*@C
885:    DSRegister - Adds a direct solver to the DS package.

887:    Not collective

889:    Input Parameters:
890: +  name - name of a new user-defined DS
891: -  routine_create - routine to create context

893:    Notes:
894:    DSRegister() may be called multiple times to add several user-defined
895:    direct solvers.

897:    Level: advanced

899: .seealso: DSRegisterAll()
900: @*/
901: PetscErrorCode DSRegister(const char *name,PetscErrorCode (*function)(DS))
902: {

906:   PetscFunctionListAdd(&DSList,name,function);
907:   return(0);
908: }

910: PETSC_EXTERN PetscErrorCode DSCreate_HEP(DS);
911: PETSC_EXTERN PetscErrorCode DSCreate_NHEP(DS);
912: PETSC_EXTERN PetscErrorCode DSCreate_GHEP(DS);
913: PETSC_EXTERN PetscErrorCode DSCreate_GHIEP(DS);
914: PETSC_EXTERN PetscErrorCode DSCreate_GNHEP(DS);
915: PETSC_EXTERN PetscErrorCode DSCreate_SVD(DS);
916: PETSC_EXTERN PetscErrorCode DSCreate_PEP(DS);
917: PETSC_EXTERN PetscErrorCode DSCreate_NEP(DS);

921: /*@C
922:    DSRegisterAll - Registers all of the direct solvers in the DS package.

924:    Not Collective

926:    Level: advanced
927: @*/
928: PetscErrorCode DSRegisterAll(void)
929: {

933:   if (DSRegisterAllCalled) return(0);
934:   DSRegisterAllCalled = PETSC_TRUE;
935:   DSRegister(DSHEP,DSCreate_HEP);
936:   DSRegister(DSNHEP,DSCreate_NHEP);
937:   DSRegister(DSGHEP,DSCreate_GHEP);
938:   DSRegister(DSGHIEP,DSCreate_GHIEP);
939:   DSRegister(DSGNHEP,DSCreate_GNHEP);
940:   DSRegister(DSSVD,DSCreate_SVD);
941:   DSRegister(DSPEP,DSCreate_PEP);
942:   DSRegister(DSNEP,DSCreate_NEP);
943:   return(0);
944: }