Actual source code: dsbasic.c

  1: /*
  2:    Basic DS routines

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

  8:    This file is part of SLEPc.
  9:       
 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: PetscFList       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"};

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

 39:    Level: developer

 41: .seealso: SlepcFinalize()
 42: @*/
 43: PetscErrorCode DSFinalizePackage(void)
 44: {
 46:   DSPackageInitialized = PETSC_FALSE;
 47:   DSList               = 0;
 48:   DSRegisterAllCalled  = PETSC_FALSE;
 49:   return(0);
 50: }

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

 59:   Input Parameter:
 60:   path - The dynamic library path, or PETSC_NULL

 62:   Level: developer

 64: .seealso: SlepcInitialize()
 65: @*/
 66: PetscErrorCode DSInitializePackage(const char *path)
 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(path);
 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(PETSC_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(PETSC_NULL,"-log_summary_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: /*@C
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:   SlepcHeaderCreate(ds,_p_DS,struct _DSOps,DS_CLASSID,-1,"DS","Direct Solver (or Dense System)","DS",comm,DSDestroy,DSView);
134:   *newds       = ds;
135:   ds->state    = DS_STATE_RAW;
136:   ds->method   = 0;
137:   ds->compact  = PETSC_FALSE;
138:   ds->refined  = PETSC_FALSE;
139:   ds->extrarow = PETSC_FALSE;
140:   ds->ld       = 0;
141:   ds->l        = 0;
142:   ds->n        = 0;
143:   ds->m        = 0;
144:   ds->k        = 0;
145:   ds->t        = 0;
146:   for (i=0;i<DS_NUM_MAT;i++) {
147:     ds->mat[i]  = PETSC_NULL;
148:     ds->rmat[i] = PETSC_NULL;
149:   }
150:   ds->perm     = PETSC_NULL;
151:   ds->work     = PETSC_NULL;
152:   ds->rwork    = PETSC_NULL;
153:   ds->iwork    = PETSC_NULL;
154:   ds->lwork    = 0;
155:   ds->lrwork   = 0;
156:   ds->liwork   = 0;
157:   ds->comp_fun = PETSC_NULL;
158:   ds->comp_ctx = PETSC_NULL;
159:   return(0);
160: }

164: /*@C
165:    DSSetOptionsPrefix - Sets the prefix used for searching for all 
166:    DS options in the database.

168:    Logically Collective on DS

170:    Input Parameters:
171: +  ds - the direct solver context
172: -  prefix - the prefix string to prepend to all DS option requests

174:    Notes:
175:    A hyphen (-) must NOT be given at the beginning of the prefix name.
176:    The first character of all runtime options is AUTOMATICALLY the
177:    hyphen.

179:    Level: advanced

181: .seealso: DSAppendOptionsPrefix()
182: @*/
183: PetscErrorCode DSSetOptionsPrefix(DS ds,const char *prefix)
184: {

189:   PetscObjectSetOptionsPrefix((PetscObject)ds,prefix);
190:   return(0);
191: }

195: /*@C
196:    DSAppendOptionsPrefix - Appends to the prefix used for searching for all 
197:    DS options in the database.

199:    Logically Collective on DS

201:    Input Parameters:
202: +  ds - the direct solver context
203: -  prefix - the prefix string to prepend to all DS option requests

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

209:    Level: advanced

211: .seealso: DSSetOptionsPrefix()
212: @*/
213: PetscErrorCode DSAppendOptionsPrefix(DS ds,const char *prefix)
214: {

219:   PetscObjectAppendOptionsPrefix((PetscObject)ds,prefix);
220:   return(0);
221: }

225: /*@C
226:    DSGetOptionsPrefix - Gets the prefix used for searching for all 
227:    DS options in the database.

229:    Not Collective

231:    Input Parameters:
232: .  ds - the direct solver context

234:    Output Parameters:
235: .  prefix - pointer to the prefix string used is returned

237:    Notes: On the fortran side, the user should pass in a string 'prefix' of
238:    sufficient length to hold the prefix.

240:    Level: advanced

242: .seealso: DSSetOptionsPrefix(), DSAppendOptionsPrefix()
243: @*/
244: PetscErrorCode DSGetOptionsPrefix(DS ds,const char *prefix[])
245: {

251:   PetscObjectGetOptionsPrefix((PetscObject)ds,prefix);
252:   return(0);
253: }

257: /*@C
258:    DSSetType - Selects the type for the DS object.

260:    Logically Collective on DS

262:    Input Parameter:
263: +  ds   - the direct solver context
264: -  type - a known type

266:    Level: intermediate

268: .seealso: DSGetType()
269: @*/
270: PetscErrorCode DSSetType(DS ds,const DSType type)
271: {
272:   PetscErrorCode ierr,(*r)(DS);
273:   PetscBool      match;


279:   PetscObjectTypeCompare((PetscObject)ds,type,&match);
280:   if (match) return(0);

282:    PetscFListFind(DSList,((PetscObject)ds)->comm,type,PETSC_TRUE,(void (**)(void))&r);
283:   if (!r) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested DS type %s",type);

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

287:   PetscObjectChangeTypeName((PetscObject)ds,type);
288:   (*r)(ds);
289:   return(0);
290: }

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

297:    Not Collective

299:    Input Parameter:
300: .  ds - the direct solver context

302:    Output Parameter:
303: .  name - name of the direct solver

305:    Level: intermediate

307: .seealso: DSSetType()
308: @*/
309: PetscErrorCode DSGetType(DS ds,const DSType *type)
310: {
314:   *type = ((PetscObject)ds)->type_name;
315:   return(0);
316: }

320: /*@
321:    DSSetMethod - Selects the method to be used to solve the problem.

323:    Logically Collective on DS

325:    Input Parameter:
326: +  ds   - the direct solver context
327: -  meth - an index indentifying the method

329:    Level: intermediate

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

346: /*@
347:    DSGetMethod - Gets the method currently used in the DS.

349:    Not Collective

351:    Input Parameter:
352: .  ds - the direct solver context

354:    Output Parameter:
355: .  meth - identifier of the method

357:    Level: intermediate

359: .seealso: DSSetMethod()
360: @*/
361: PetscErrorCode DSGetMethod(DS ds,PetscInt *meth)
362: {
366:   *meth = ds->method;
367:   return(0);
368: }

372: /*@
373:    DSSetCompact - Switch to compact storage of matrices.

375:    Logically Collective on DS

377:    Input Parameter:
378: +  ds   - the direct solver context
379: -  comp - a boolean flag

381:    Notes:
382:    Compact storage is used in some DS types such as DSHEP when the matrix
383:    is tridiagonal. This flag can be used to indicate whether the user
384:    provides the matrix entries via the compact form (the tridiagonal DS_MAT_T)
385:    or the non-compact one (DS_MAT_A).

387:    The default is PETSC_FALSE.

389:    Level: advanced

391: .seealso: DSGetCompact()
392: @*/
393: PetscErrorCode DSSetCompact(DS ds,PetscBool comp)
394: {
398:   ds->compact = comp;
399:   return(0);
400: }

404: /*@
405:    DSGetCompact - Gets the compact storage flag.

407:    Not Collective

409:    Input Parameter:
410: .  ds - the direct solver context

412:    Output Parameter:
413: .  comp - the flag

415:    Level: advanced

417: .seealso: DSSetCompact()
418: @*/
419: PetscErrorCode DSGetCompact(DS ds,PetscBool *comp)
420: {
424:   *comp = ds->compact;
425:   return(0);
426: }

430: /*@
431:    DSSetExtraRow - Sets a flag to indicate that the matrix has one extra
432:    row.

434:    Logically Collective on DS

436:    Input Parameter:
437: +  ds  - the direct solver context
438: -  ext - a boolean flag

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

446:    The default is PETSC_FALSE.

448:    Level: advanced

450: .seealso: DSSolve(), DSAllocate(), DSGetExtraRow()
451: @*/
452: PetscErrorCode DSSetExtraRow(DS ds,PetscBool ext)
453: {
457:   if (ds->n>0 && ds->n==ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Cannot set extra row after setting n=ld");
458:   ds->extrarow = ext;
459:   return(0);
460: }

464: /*@
465:    DSGetExtraRow - Gets the extra row flag.

467:    Not Collective

469:    Input Parameter:
470: .  ds - the direct solver context

472:    Output Parameter:
473: .  ext - the flag

475:    Level: advanced

477: .seealso: DSSetExtraRow()
478: @*/
479: PetscErrorCode DSGetExtraRow(DS ds,PetscBool *ext)
480: {
484:   *ext = ds->extrarow;
485:   return(0);
486: }

490: /*@
491:    DSSetRefined - Sets a flag to indicate that refined vectors must be
492:    computed.

494:    Logically Collective on DS

496:    Input Parameter:
497: +  ds  - the direct solver context
498: -  ref - a boolean flag

500:    Notes:
501:    Normally the vectors returned in DS_MAT_X are eigenvectors of the
502:    projected matrix. With this flag activated, DSVectors() will return
503:    the right singular vector of the smallest singular value of matrix
504:    \tilde{A}-theta*I, where \tilde{A} is the extended (n+1)xn matrix
505:    and theta is the Ritz value. This is used in the refined Ritz
506:    approximation.

508:    The default is PETSC_FALSE.

510:    Level: advanced

512: .seealso: DSVectors(), DSGetRefined()
513: @*/
514: PetscErrorCode DSSetRefined(DS ds,PetscBool ref)
515: {
519:   ds->refined = ref;
520:   return(0);
521: }

525: /*@
526:    DSGetRefined - Gets the refined vectors flag.

528:    Not Collective

530:    Input Parameter:
531: .  ds - the direct solver context

533:    Output Parameter:
534: .  ref - the flag

536:    Level: advanced

538: .seealso: DSSetRefined()
539: @*/
540: PetscErrorCode DSGetRefined(DS ds,PetscBool *ref)
541: {
545:   *ref = ds->refined;
546:   return(0);
547: }

551: /*@C
552:    DSSetEigenvalueComparison - Specifies the eigenvalue comparison function
553:    to be used for sorting.

555:    Logically Collective on DS

557:    Input Parameters:
558: +  ds  - the direct solver context
559: .  fun - a pointer to the comparison function
560: -  ctx - a context pointer (the last parameter to the comparison function)

562:    Calling Sequence of fun:
563: $  func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)

565: +   ar     - real part of the 1st eigenvalue
566: .   ai     - imaginary part of the 1st eigenvalue
567: .   br     - real part of the 2nd eigenvalue
568: .   bi     - imaginary part of the 2nd eigenvalue
569: .   res    - result of comparison
570: -   ctx    - optional context, as set by DSSetEigenvalueComparison()

572:    Note:
573:    The returning parameter 'res' can be:
574: +  negative - if the 1st eigenvalue is preferred to the 2st one
575: .  zero     - if both eigenvalues are equally preferred
576: -  positive - if the 2st eigenvalue is preferred to the 1st one

578:    Level: developer

580: .seealso: DSSort()
581: @*/
582: PetscErrorCode DSSetEigenvalueComparison(DS ds,PetscErrorCode (*fun)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)
583: {
586:   ds->comp_fun = fun;
587:   ds->comp_ctx = ctx;
588:   return(0);
589: }

593: /*@C
594:    DSGetEigenvalueComparison - Gets the eigenvalue comparison function
595:    used for sorting.

597:    Not Collective

599:    Input Parameter:
600: .  ds  - the direct solver context

602:    Output Parameters:
603: +  fun - a pointer to the comparison function
604: -  ctx - a context pointer (the last parameter to the comparison function)

606:    Calling Sequence of fun:
607: $  func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)

609: +   ar     - real part of the 1st eigenvalue
610: .   ai     - imaginary part of the 1st eigenvalue
611: .   br     - real part of the 2nd eigenvalue
612: .   bi     - imaginary part of the 2nd eigenvalue
613: .   res    - result of comparison
614: -   ctx    - optional context, as set by DSSetEigenvalueComparison()

616:    Note:
617:    The returning parameter 'res' can be:
618: +  negative - if the 1st eigenvalue is preferred to the 2st one
619: .  zero     - if both eigenvalues are equally preferred
620: -  positive - if the 2st eigenvalue is preferred to the 1st one

622:    Level: developer

624: .seealso: DSSort(), DSSetEigenvalueComparison()
625: @*/
626: PetscErrorCode DSGetEigenvalueComparison(DS ds,PetscErrorCode (**fun)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void** ctx)
627: {
630:   if (fun) *fun = ds->comp_fun;
631:   if (ctx) *ctx = ds->comp_ctx;
632:   return(0);
633: }

637: /*@
638:    DSSetFromOptions - Sets DS options from the options database.

640:    Collective on DS

642:    Input Parameters:
643: .  ds - the direct solver context

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

648:    Level: beginner
649: @*/
650: PetscErrorCode DSSetFromOptions(DS ds)
651: {
653:   PetscInt       meth;

657:   if (!DSRegisterAllCalled) { DSRegisterAll(PETSC_NULL); }
658:   /* Set default type (we do not allow changing it with -ds_type) */
659:   if (!((PetscObject)ds)->type_name) {
660:     DSSetType(ds,DSNHEP);
661:   }
662:   PetscOptionsBegin(((PetscObject)ds)->comm,((PetscObject)ds)->prefix,"Direct Solver (DS) Options","DS");
663:     meth = 0;
664:     PetscOptionsInt("-ds_method","Method to be used for the dense system","DSSetMethod",ds->method,&meth,PETSC_NULL);
665:     DSSetMethod(ds,meth);
666:     PetscObjectProcessOptionsHandlers((PetscObject)ds);
667:   PetscOptionsEnd();
668:   return(0);
669: }

673: /*@C
674:    DSView - Prints the DS data structure.

676:    Collective on DS

678:    Input Parameters:
679: +  ds - the direct solver context
680: -  viewer - optional visualization context

682:    Note:
683:    The available visualization contexts include
684: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
685: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
686:          output where only the first processor opens
687:          the file.  All other processors send their 
688:          data to the first processor to print. 

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

693:    Level: beginner

695: .seealso: PetscViewerASCIIOpen()
696: @*/
697: PetscErrorCode DSView(DS ds,PetscViewer viewer)
698: {
699:   PetscBool         isascii,issvd;
700:   const char        *state;
701:   PetscViewerFormat format;
702:   PetscErrorCode    ierr;

706:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(((PetscObject)ds)->comm);
709:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
710:   if (isascii) {
711:     PetscViewerGetFormat(viewer,&format);
712:     PetscObjectPrintClassNamePrefixType((PetscObject)ds,viewer,"DS Object");
713:     if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
714:       switch (ds->state) {
715:         case DS_STATE_RAW:          state = "raw"; break;
716:         case DS_STATE_INTERMEDIATE: state = "intermediate"; break;
717:         case DS_STATE_CONDENSED:    state = "condensed"; break;
718:         case DS_STATE_TRUNCATED:    state = "truncated"; break;
719:         default: SETERRQ(((PetscObject)ds)->comm,1,"Wrong value of ds->state");
720:       }
721:       PetscViewerASCIIPrintf(viewer,"  current state: %s\n",state);
722:       PetscObjectTypeCompare((PetscObject)ds,DSSVD,&issvd);
723:       if (issvd) {
724:         PetscViewerASCIIPrintf(viewer,"  dimensions: ld=%d, n=%d, m=%d, l=%d, k=%d",ds->ld,ds->n,ds->m,ds->l,ds->k);
725:       } else {
726:         PetscViewerASCIIPrintf(viewer,"  dimensions: ld=%d, n=%d, l=%d, k=%d",ds->ld,ds->n,ds->l,ds->k);
727:       }
728:       if (ds->state==DS_STATE_TRUNCATED) {
729:         PetscViewerASCIIPrintf(viewer,", t=%d\n",ds->t);
730:       } else {
731:         PetscViewerASCIIPrintf(viewer,"\n");
732:       }
733:       PetscViewerASCIIPrintf(viewer,"  flags:%s%s%s\n",ds->compact?" compact":"",ds->extrarow?" extrarow":"",ds->refined?" refined":"");
734:     }
735:     if (ds->ops->view) {
736:       PetscViewerASCIIPushTab(viewer);
737:       (*ds->ops->view)(ds,viewer);
738:       PetscViewerASCIIPopTab(viewer);
739:     }
740:   } else SETERRQ1(((PetscObject)ds)->comm,1,"Viewer type %s not supported for DS",((PetscObject)viewer)->type_name);
741:   return(0);
742: }

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

749:    Logically Collective on DS

751:    Input Parameters:
752: +  ds - the direct solver context
753: -  ld - leading dimension (maximum allowed dimension for the matrices, including
754:         the extra row if present)

756:    Level: intermediate

758: .seealso: DSGetLeadingDimension(), DSSetDimensions(), DSSetExtraRow()
759: @*/
760: PetscErrorCode DSAllocate(DS ds,PetscInt ld)
761: {

767:   if (ld<1) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Leading dimension should be at least one");
768:   ds->ld = ld;
769:   (*ds->ops->allocate)(ds,ld);
770:   return(0);
771: }

775: /*@
776:    DSReset - Resets the DS context to the initial state.

778:    Collective on DS

780:    Input Parameter:
781: .  ds - the direct solver context

783:    Level: advanced

785: .seealso: DSDestroy()
786: @*/
787: PetscErrorCode DSReset(DS ds)
788: {
789:   PetscInt       i;

794:   ds->state    = DS_STATE_RAW;
795:   ds->compact  = PETSC_FALSE;
796:   ds->refined  = PETSC_FALSE;
797:   ds->extrarow = PETSC_FALSE;
798:   ds->ld       = 0;
799:   ds->l        = 0;
800:   ds->n        = 0;
801:   ds->m        = 0;
802:   ds->k        = 0;
803:   for (i=0;i<DS_NUM_MAT;i++) {
804:     PetscFree(ds->mat[i]);
805:     PetscFree(ds->rmat[i]);
806:   }
807:   PetscFree(ds->perm);
808:   PetscFree(ds->work);
809:   PetscFree(ds->rwork);
810:   PetscFree(ds->iwork);
811:   ds->lwork    = 0;
812:   ds->lrwork   = 0;
813:   ds->liwork   = 0;
814:   ds->comp_fun = PETSC_NULL;
815:   ds->comp_ctx = PETSC_NULL;
816:   return(0);
817: }

821: /*@C
822:    DSDestroy - Destroys DS context that was created with DSCreate().

824:    Collective on DS

826:    Input Parameter:
827: .  ds - the direct solver context

829:    Level: beginner

831: .seealso: DSCreate()
832: @*/
833: PetscErrorCode DSDestroy(DS *ds)
834: {

838:   if (!*ds) return(0);
840:   if (--((PetscObject)(*ds))->refct > 0) { *ds = 0; return(0); }
841:   DSReset(*ds);
842:   PetscHeaderDestroy(ds);
843:   return(0);
844: }

848: /*@C
849:    DSRegister - See DSRegisterDynamic()

851:    Level: advanced
852: @*/
853: PetscErrorCode DSRegister(const char *sname,const char *path,const char *name,PetscErrorCode (*function)(DS))
854: {
856:   char           fullname[PETSC_MAX_PATH_LEN];

859:   PetscFListConcat(path,name,fullname);
860:   PetscFListAdd(&DSList,sname,fullname,(void (*)(void))function);
861:   return(0);
862: }

866: /*@
867:    DSRegisterDestroy - Frees the list of DS methods that were
868:    registered by DSRegisterDynamic().

870:    Not Collective

872:    Level: advanced

874: .seealso: DSRegisterDynamic(), DSRegisterAll()
875: @*/
876: PetscErrorCode DSRegisterDestroy(void)
877: {

881:   PetscFListDestroy(&DSList);
882:   DSRegisterAllCalled = PETSC_FALSE;
883:   return(0);
884: }

886: EXTERN_C_BEGIN
887: extern PetscErrorCode DSCreate_HEP(DS);
888: extern PetscErrorCode DSCreate_NHEP(DS);
889: extern PetscErrorCode DSCreate_GHEP(DS);
890: extern PetscErrorCode DSCreate_GHIEP(DS);
891: extern PetscErrorCode DSCreate_GNHEP(DS);
892: extern PetscErrorCode DSCreate_SVD(DS);
893: EXTERN_C_END

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

900:    Not Collective

902:    Input Parameter:
903: .  path - the library where the routines are to be found (optional)

905:    Level: advanced
906: @*/
907: PetscErrorCode DSRegisterAll(const char *path)
908: {

912:   DSRegisterAllCalled = PETSC_TRUE;
913:   DSRegisterDynamic(DSHEP,path,"DSCreate_HEP",DSCreate_HEP);
914:   DSRegisterDynamic(DSNHEP,path,"DSCreate_NHEP",DSCreate_NHEP);
915:   DSRegisterDynamic(DSGHEP,path,"DSCreate_GHEP",DSCreate_GHEP);
916:   DSRegisterDynamic(DSGHIEP,path,"DSCreate_GHIEP",DSCreate_GHIEP);
917:   DSRegisterDynamic(DSGNHEP,path,"DSCreate_GNHEP",DSCreate_GNHEP);
918:   DSRegisterDynamic(DSSVD,path,"DSCreate_SVD",DSCreate_SVD);
919:   return(0);
920: }