Actual source code: dsbasic.c

slepc-3.5.2 2014-10-10
Report Typos and Errors
  1: /*
  2:    Basic DS routines

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2014, 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_Function = 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","F","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("DSFunction",DS_CLASSID,&DS_Function);
 83:   PetscLogEventRegister("DSVectors",DS_CLASSID,&DS_Vectors);
 84:   PetscLogEventRegister("DSOther",DS_CLASSID,&DS_Other);
 85:   /* Process info exclusions */
 86:   PetscOptionsGetString(NULL,"-info_exclude",logList,256,&opt);
 87:   if (opt) {
 88:     PetscStrstr(logList,"ds",&className);
 89:     if (className) {
 90:       PetscInfoDeactivateClass(DS_CLASSID);
 91:     }
 92:   }
 93:   /* Process summary exclusions */
 94:   PetscOptionsGetString(NULL,"-log_summary_exclude",logList,256,&opt);
 95:   if (opt) {
 96:     PetscStrstr(logList,"ds",&className);
 97:     if (className) {
 98:       PetscLogEventDeactivateClass(DS_CLASSID);
 99:     }
100:   }
101:   PetscRegisterFinalize(DSFinalizePackage);
102:   return(0);
103: }

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

110:    Collective on MPI_Comm

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

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

118:    Level: beginner

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

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

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

138:   ds->state         = DS_STATE_RAW;
139:   ds->method        = 0;
140:   ds->funmethod     = 0;
141:   ds->compact       = PETSC_FALSE;
142:   ds->refined       = PETSC_FALSE;
143:   ds->extrarow      = PETSC_FALSE;
144:   ds->ld            = 0;
145:   ds->l             = 0;
146:   ds->n             = 0;
147:   ds->m             = 0;
148:   ds->k             = 0;
149:   ds->t             = 0;
150:   ds->bs            = 1;
151:   ds->nf            = 0;
152:   for (i=0;i<DS_NUM_EXTRA;i++) ds->f[i] = NULL;
153:   ds->sc            = NULL;

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

168:   *newds = ds;
169:   return(0);
170: }

174: /*@C
175:    DSSetOptionsPrefix - Sets the prefix used for searching for all
176:    DS options in the database.

178:    Logically Collective on DS

180:    Input Parameters:
181: +  ds - the direct solver context
182: -  prefix - the prefix string to prepend to all DS option requests

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

189:    Level: advanced

191: .seealso: DSAppendOptionsPrefix()
192: @*/
193: PetscErrorCode DSSetOptionsPrefix(DS ds,const char *prefix)
194: {

199:   PetscObjectSetOptionsPrefix((PetscObject)ds,prefix);
200:   return(0);
201: }

205: /*@C
206:    DSAppendOptionsPrefix - Appends to the prefix used for searching for all
207:    DS options in the database.

209:    Logically Collective on DS

211:    Input Parameters:
212: +  ds - the direct solver context
213: -  prefix - the prefix string to prepend to all DS option requests

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

219:    Level: advanced

221: .seealso: DSSetOptionsPrefix()
222: @*/
223: PetscErrorCode DSAppendOptionsPrefix(DS ds,const char *prefix)
224: {

229:   PetscObjectAppendOptionsPrefix((PetscObject)ds,prefix);
230:   return(0);
231: }

235: /*@C
236:    DSGetOptionsPrefix - Gets the prefix used for searching for all
237:    DS options in the database.

239:    Not Collective

241:    Input Parameters:
242: .  ds - the direct solver context

244:    Output Parameters:
245: .  prefix - pointer to the prefix string used is returned

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

250:    Level: advanced

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

261:   PetscObjectGetOptionsPrefix((PetscObject)ds,prefix);
262:   return(0);
263: }

267: /*@C
268:    DSSetType - Selects the type for the DS object.

270:    Logically Collective on DS

272:    Input Parameter:
273: +  ds   - the direct solver context
274: -  type - a known type

276:    Level: intermediate

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


289:   PetscObjectTypeCompare((PetscObject)ds,type,&match);
290:   if (match) return(0);

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

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

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

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

307:    Not Collective

309:    Input Parameter:
310: .  ds - the direct solver context

312:    Output Parameter:
313: .  name - name of the direct solver

315:    Level: intermediate

317: .seealso: DSSetType()
318: @*/
319: PetscErrorCode DSGetType(DS ds,DSType *type)
320: {
324:   *type = ((PetscObject)ds)->type_name;
325:   return(0);
326: }

330: /*@
331:    DSSetMethod - Selects the method to be used to solve the problem.

333:    Logically Collective on DS

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

339:    Level: intermediate

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

356: /*@
357:    DSGetMethod - Gets the method currently used in the DS.

359:    Not Collective

361:    Input Parameter:
362: .  ds - the direct solver context

364:    Output Parameter:
365: .  meth - identifier of the method

367:    Level: intermediate

369: .seealso: DSSetMethod()
370: @*/
371: PetscErrorCode DSGetMethod(DS ds,PetscInt *meth)
372: {
376:   *meth = ds->method;
377:   return(0);
378: }

382: /*@
383:    DSSetFunctionMethod - Selects the method to be used to compute a matrix function.

385:    Logically Collective on DS

387:    Input Parameter:
388: +  ds   - the direct solver context
389: -  meth - an index indentifying the function method

391:    Level: intermediate

393: .seealso: DSGetFunctionMethod()
394: @*/
395: PetscErrorCode DSSetFunctionMethod(DS ds,PetscInt meth)
396: {
400:   if (meth<0) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"The method must be a non-negative integer");
401:   if (meth>DS_MAX_FUN) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"Too large value for the method");
402:   ds->funmethod = meth;
403:   return(0);
404: }

408: /*@
409:    DSGetFunctionMethod - Gets the method currently used to compute a matrix function.

411:    Not Collective

413:    Input Parameter:
414: .  ds - the direct solver context

416:    Output Parameter:
417: .  meth - identifier of the function method

419:    Level: intermediate

421: .seealso: DSSetFunctionMethod()
422: @*/
423: PetscErrorCode DSGetFunctionMethod(DS ds,PetscInt *meth)
424: {
428:   *meth = ds->funmethod;
429:   return(0);
430: }

434: /*@
435:    DSSetCompact - Switch to compact storage of matrices.

437:    Logically Collective on DS

439:    Input Parameter:
440: +  ds   - the direct solver context
441: -  comp - a boolean flag

443:    Notes:
444:    Compact storage is used in some DS types such as DSHEP when the matrix
445:    is tridiagonal. This flag can be used to indicate whether the user
446:    provides the matrix entries via the compact form (the tridiagonal DS_MAT_T)
447:    or the non-compact one (DS_MAT_A).

449:    The default is PETSC_FALSE.

451:    Level: advanced

453: .seealso: DSGetCompact()
454: @*/
455: PetscErrorCode DSSetCompact(DS ds,PetscBool comp)
456: {
460:   ds->compact = comp;
461:   return(0);
462: }

466: /*@
467:    DSGetCompact - Gets the compact storage flag.

469:    Not Collective

471:    Input Parameter:
472: .  ds - the direct solver context

474:    Output Parameter:
475: .  comp - the flag

477:    Level: advanced

479: .seealso: DSSetCompact()
480: @*/
481: PetscErrorCode DSGetCompact(DS ds,PetscBool *comp)
482: {
486:   *comp = ds->compact;
487:   return(0);
488: }

492: /*@
493:    DSSetExtraRow - Sets a flag to indicate that the matrix has one extra
494:    row.

496:    Logically Collective on DS

498:    Input Parameter:
499: +  ds  - the direct solver context
500: -  ext - a boolean flag

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

508:    The default is PETSC_FALSE.

510:    Level: advanced

512: .seealso: DSSolve(), DSAllocate(), DSGetExtraRow()
513: @*/
514: PetscErrorCode DSSetExtraRow(DS ds,PetscBool ext)
515: {
519:   if (ds->n>0 && ds->n==ds->ld) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ORDER,"Cannot set extra row after setting n=ld");
520:   ds->extrarow = ext;
521:   return(0);
522: }

526: /*@
527:    DSGetExtraRow - Gets the extra row flag.

529:    Not Collective

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

534:    Output Parameter:
535: .  ext - the flag

537:    Level: advanced

539: .seealso: DSSetExtraRow()
540: @*/
541: PetscErrorCode DSGetExtraRow(DS ds,PetscBool *ext)
542: {
546:   *ext = ds->extrarow;
547:   return(0);
548: }

552: /*@
553:    DSSetRefined - Sets a flag to indicate that refined vectors must be
554:    computed.

556:    Logically Collective on DS

558:    Input Parameter:
559: +  ds  - the direct solver context
560: -  ref - a boolean flag

562:    Notes:
563:    Normally the vectors returned in DS_MAT_X are eigenvectors of the
564:    projected matrix. With this flag activated, DSVectors() will return
565:    the right singular vector of the smallest singular value of matrix
566:    \tilde{A}-theta*I, where \tilde{A} is the extended (n+1)xn matrix
567:    and theta is the Ritz value. This is used in the refined Ritz
568:    approximation.

570:    The default is PETSC_FALSE.

572:    Level: advanced

574: .seealso: DSVectors(), DSGetRefined()
575: @*/
576: PetscErrorCode DSSetRefined(DS ds,PetscBool ref)
577: {
581:   ds->refined = ref;
582:   return(0);
583: }

587: /*@
588:    DSGetRefined - Gets the refined vectors flag.

590:    Not Collective

592:    Input Parameter:
593: .  ds - the direct solver context

595:    Output Parameter:
596: .  ref - the flag

598:    Level: advanced

600: .seealso: DSSetRefined()
601: @*/
602: PetscErrorCode DSGetRefined(DS ds,PetscBool *ref)
603: {
607:   *ref = ds->refined;
608:   return(0);
609: }

613: /*@
614:    DSSetBlockSize - Sets the block size.

616:    Logically Collective on DS

618:    Input Parameter:
619: +  ds - the direct solver context
620: -  bs - the block size

622:    Level: intermediate

624: .seealso: DSGetBlockSize()
625: @*/
626: PetscErrorCode DSSetBlockSize(DS ds,PetscInt bs)
627: {
631:   if (bs<1) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"The block size must be at least one");
632:   ds->bs = bs;
633:   return(0);
634: }

638: /*@
639:    DSGetBlockSize - Gets the block size.

641:    Not Collective

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

646:    Output Parameter:
647: .  bs - block size

649:    Level: intermediate

651: .seealso: DSSetBlockSize()
652: @*/
653: PetscErrorCode DSGetBlockSize(DS ds,PetscInt *bs)
654: {
658:   *bs = ds->bs;
659:   return(0);
660: }

664: /*@C
665:    DSSetSlepcSC - Sets the sorting criterion context.

667:    Not Collective

669:    Input Parameters:
670: +  ds - the direct solver context
671: -  sc - a pointer to the sorting criterion context

673:    Level: developer

675: .seealso: DSGetSlepcSC(), DSSort()
676: @*/
677: PetscErrorCode DSSetSlepcSC(DS ds,SlepcSC sc)
678: {

684:   if (ds->sc) {
685:     PetscFree(ds->sc);
686:   }
687:   ds->sc = sc;
688:   return(0);
689: }

693: /*@C
694:    DSGetSlepcSC - Gets the sorting criterion context.

696:    Not Collective

698:    Input Parameter:
699: .  ds - the direct solver context

701:    Output Parameters:
702: .  sc - a pointer to the sorting criterion context

704:    Level: developer

706: .seealso: DSSetSlepcSC(), DSSort()
707: @*/
708: PetscErrorCode DSGetSlepcSC(DS ds,SlepcSC *sc)
709: {

715:   if (!ds->sc) {
716:     PetscNewLog(ds,&ds->sc);
717:   }
718:   *sc = ds->sc;
719:   return(0);
720: }

724: /*@
725:    DSSetFN - Sets a number of functions to be used internally by DS.

727:    Collective on DS and FN

729:    Input Parameters:
730: +  ds - the direct solver context
731: .  n  - number of functions
732: -  f  - array of functions

734:    Notes:
735:    In the basic usage, only one function is used, for instance to
736:    evaluate a function of the projected matrix. In the context of nonlinear
737:    eigensolvers, there are as many functions as terms in the split
738:    nonlinear operator T(lambda) = sum_i A_i*f_i(lambda).

740:    This function must be called before DSAllocate(). Then DSAllocate()
741:    will allocate an extra matrix per each function.

743:    Level: developer

745: .seealso: DSGetFN(), DSGetFN(), DSAllocate()
746:  @*/
747: PetscErrorCode DSSetFN(DS ds,PetscInt n,FN f[])
748: {
749:   PetscInt       i;

755:   if (n<=0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more functions, you have %D",n);
756:   if (n>DS_NUM_EXTRA) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Too many functions, you specified %D but the limit is %D",n,DS_NUM_EXTRA);
757:   if (ds->ld) { PetscInfo(ds,"DSSetFN() called after DSAllocate()\n"); }
760:   for (i=0;i<ds->nf;i++) {
761:     FNDestroy(&ds->f[i]);
762:   }
763:   for (i=0;i<n;i++) {
765:     PetscObjectReference((PetscObject)f[i]);
766:     ds->f[i] = f[i];
767:   }
768:   ds->nf = n;
769:   return(0);
770: }

774: /*@
775:    DSGetFN - Gets the functions associated with this DS.

777:    Not collective, though parallel FNs are returned if the DS is parallel

779:    Input Parameter:
780: +  ds - the direct olver context
781: -  k  - the index of the requested function (starting in 0)

783:    Output Parameter:
784: .  f - the function

786:    Level: developer

788: .seealso: DSSetFN()
789: @*/
790: PetscErrorCode DSGetFN(DS ds,PetscInt k,FN *f)
791: {
794:   if (k<0 || k>=ds->nf) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %d",ds->nf-1);
796:   *f = ds->f[k];
797:   return(0);
798: }

802: /*@
803:    DSGetNumFN - Returns the number of functions stored internally by
804:    the DS.

806:    Not collective

808:    Input Parameter:
809: .  ds - the direct solver context

811:    Output Parameters:
812: .  n - the number of functions passed in DSSetFN()

814:    Level: developer

816: .seealso: DSSetFN()
817: @*/
818: PetscErrorCode DSGetNumFN(DS ds,PetscInt *n)
819: {
823:   *n = ds->nf;
824:   return(0);
825: }

829: /*@
830:    DSSetFromOptions - Sets DS options from the options database.

832:    Collective on DS

834:    Input Parameters:
835: .  ds - the direct solver context

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

840:    Level: beginner
841: @*/
842: PetscErrorCode DSSetFromOptions(DS ds)
843: {
845:   PetscInt       bs,meth;
846:   PetscBool      flag;

850:   if (!DSRegisterAllCalled) { DSRegisterAll(); }
851:   /* Set default type (we do not allow changing it with -ds_type) */
852:   if (!((PetscObject)ds)->type_name) {
853:     DSSetType(ds,DSNHEP);
854:   }
855:   PetscObjectOptionsBegin((PetscObject)ds);
856:     PetscOptionsInt("-ds_block_size","Block size for the dense system solver","DSSetBlockSize",ds->bs,&bs,&flag);
857:     if (flag) { DSSetBlockSize(ds,bs); }
858:     PetscOptionsInt("-ds_method","Method to be used for the dense system","DSSetMethod",ds->method,&meth,&flag);
859:     if (flag) { DSSetMethod(ds,meth); }
860:     PetscOptionsInt("-ds_function_method","Method to be used to compute a matrix function","DSSetFunctionMethod",ds->funmethod,&meth,&flag);
861:     if (flag) { DSSetFunctionMethod(ds,meth); }
862:     PetscObjectProcessOptionsHandlers((PetscObject)ds);
863:   PetscOptionsEnd();
864:   return(0);
865: }

869: /*@C
870:    DSView - Prints the DS data structure.

872:    Collective on DS

874:    Input Parameters:
875: +  ds - the direct solver context
876: -  viewer - optional visualization context

878:    Note:
879:    The available visualization contexts include
880: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
881: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
882:          output where only the first processor opens
883:          the file.  All other processors send their
884:          data to the first processor to print.

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

889:    Level: beginner

891: .seealso: DSViewMat()
892: @*/
893: PetscErrorCode DSView(DS ds,PetscViewer viewer)
894: {
895:   PetscBool         isascii,issvd;
896:   const char        *state;
897:   PetscViewerFormat format;
898:   PetscErrorCode    ierr;

902:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ds));
905:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
906:   if (isascii) {
907:     PetscViewerGetFormat(viewer,&format);
908:     PetscObjectPrintClassNamePrefixType((PetscObject)ds,viewer);
909:     if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
910:       switch (ds->state) {
911:         case DS_STATE_RAW:          state = "raw"; break;
912:         case DS_STATE_INTERMEDIATE: state = "intermediate"; break;
913:         case DS_STATE_CONDENSED:    state = "condensed"; break;
914:         case DS_STATE_TRUNCATED:    state = "truncated"; break;
915:         default: SETERRQ(PetscObjectComm((PetscObject)ds),1,"Wrong value of ds->state");
916:       }
917:       PetscViewerASCIIPrintf(viewer,"  current state: %s\n",state);
918:       PetscObjectTypeCompare((PetscObject)ds,DSSVD,&issvd);
919:       if (issvd) {
920:         PetscViewerASCIIPrintf(viewer,"  dimensions: ld=%D, n=%D, m=%D, l=%D, k=%D",ds->ld,ds->n,ds->m,ds->l,ds->k);
921:       } else {
922:         PetscViewerASCIIPrintf(viewer,"  dimensions: ld=%D, n=%D, l=%D, k=%D",ds->ld,ds->n,ds->l,ds->k);
923:       }
924:       if (ds->state==DS_STATE_TRUNCATED) {
925:         PetscViewerASCIIPrintf(viewer,", t=%D\n",ds->t);
926:       } else {
927:         PetscViewerASCIIPrintf(viewer,"\n");
928:       }
929:       PetscViewerASCIIPrintf(viewer,"  flags:%s%s%s\n",ds->compact?" compact":"",ds->extrarow?" extrarow":"",ds->refined?" refined":"");
930:       if (ds->nf) {
931:         PetscViewerASCIIPrintf(viewer,"  number of functions: %D\n",ds->nf);
932:       }
933:     }
934:     if (ds->ops->view) {
935:       PetscViewerASCIIPushTab(viewer);
936:       (*ds->ops->view)(ds,viewer);
937:       PetscViewerASCIIPopTab(viewer);
938:     }
939:   }
940:   return(0);
941: }

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

948:    Logically Collective on DS

950:    Input Parameters:
951: +  ds - the direct solver context
952: -  ld - leading dimension (maximum allowed dimension for the matrices, including
953:         the extra row if present)

955:    Level: intermediate

957: .seealso: DSGetLeadingDimension(), DSSetDimensions(), DSSetExtraRow()
958: @*/
959: PetscErrorCode DSAllocate(DS ds,PetscInt ld)
960: {
962:   PetscInt       i;

967:   if (ld<1) SETERRQ(PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"Leading dimension should be at least one");
968:   ds->ld = ld;
969:   (*ds->ops->allocate)(ds,ld);
970:   for (i=0;i<ds->nf;i++) {
971:     DSAllocateMat_Private(ds,DSMatExtra[i]);
972:   }
973:   return(0);
974: }

978: /*@
979:    DSReset - Resets the DS context to the initial state.

981:    Collective on DS

983:    Input Parameter:
984: .  ds - the direct solver context

986:    Level: advanced

988: .seealso: DSDestroy()
989: @*/
990: PetscErrorCode DSReset(DS ds)
991: {
992:   PetscInt       i;

997:   ds->state    = DS_STATE_RAW;
998:   ds->compact  = PETSC_FALSE;
999:   ds->refined  = PETSC_FALSE;
1000:   ds->extrarow = PETSC_FALSE;
1001:   ds->ld       = 0;
1002:   ds->l        = 0;
1003:   ds->n        = 0;
1004:   ds->m        = 0;
1005:   ds->k        = 0;
1006:   for (i=0;i<DS_NUM_MAT;i++) {
1007:     PetscFree(ds->mat[i]);
1008:     PetscFree(ds->rmat[i]);
1009:     MatDestroy(&ds->omat[i]);
1010:   }
1011:   for (i=0;i<ds->nf;i++) {
1012:     FNDestroy(&ds->f[i]);
1013:   }
1014:   ds->nf            = 0;
1015:   PetscFree(ds->perm);
1016:   PetscFree(ds->work);
1017:   PetscFree(ds->rwork);
1018:   PetscFree(ds->iwork);
1019:   ds->lwork         = 0;
1020:   ds->lrwork        = 0;
1021:   ds->liwork        = 0;
1022:   return(0);
1023: }

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

1030:    Collective on DS

1032:    Input Parameter:
1033: .  ds - the direct solver context

1035:    Level: beginner

1037: .seealso: DSCreate()
1038: @*/
1039: PetscErrorCode DSDestroy(DS *ds)
1040: {

1044:   if (!*ds) return(0);
1046:   if (--((PetscObject)(*ds))->refct > 0) { *ds = 0; return(0); }
1047:   DSReset(*ds);
1048:   PetscFree((*ds)->sc);
1049:   PetscHeaderDestroy(ds);
1050:   return(0);
1051: }

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

1058:    Not collective

1060:    Input Parameters:
1061: +  name - name of a new user-defined DS
1062: -  routine_create - routine to create context

1064:    Notes:
1065:    DSRegister() may be called multiple times to add several user-defined
1066:    direct solvers.

1068:    Level: advanced

1070: .seealso: DSRegisterAll()
1071: @*/
1072: PetscErrorCode DSRegister(const char *name,PetscErrorCode (*function)(DS))
1073: {

1077:   PetscFunctionListAdd(&DSList,name,function);
1078:   return(0);
1079: }

1081: PETSC_EXTERN PetscErrorCode DSCreate_HEP(DS);
1082: PETSC_EXTERN PetscErrorCode DSCreate_NHEP(DS);
1083: PETSC_EXTERN PetscErrorCode DSCreate_GHEP(DS);
1084: PETSC_EXTERN PetscErrorCode DSCreate_GHIEP(DS);
1085: PETSC_EXTERN PetscErrorCode DSCreate_GNHEP(DS);
1086: PETSC_EXTERN PetscErrorCode DSCreate_SVD(DS);
1087: PETSC_EXTERN PetscErrorCode DSCreate_NEP(DS);

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

1094:    Not Collective

1096:    Level: advanced
1097: @*/
1098: PetscErrorCode DSRegisterAll(void)
1099: {

1103:   DSRegisterAllCalled = PETSC_TRUE;
1104:   DSRegister(DSHEP,DSCreate_HEP);
1105:   DSRegister(DSNHEP,DSCreate_NHEP);
1106:   DSRegister(DSGHEP,DSCreate_GHEP);
1107:   DSRegister(DSGHIEP,DSCreate_GHIEP);
1108:   DSRegister(DSGNHEP,DSCreate_GNHEP);
1109:   DSRegister(DSSVD,DSCreate_SVD);
1110:   DSRegister(DSNEP,DSCreate_NEP);
1111:   return(0);
1112: }