Actual source code: dsbasic.c
slepc-3.21.0 2024-03-30
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: Basic DS routines
12: */
14: #include <slepc/private/dsimpl.h>
16: PetscFunctionList DSList = NULL;
17: PetscBool DSRegisterAllCalled = PETSC_FALSE;
18: PetscClassId DS_CLASSID = 0;
19: PetscLogEvent DS_Solve = 0,DS_Vectors = 0,DS_Synchronize = 0,DS_Other = 0;
20: static PetscBool DSPackageInitialized = PETSC_FALSE;
22: const char *DSStateTypes[] = {"RAW","INTERMEDIATE","CONDENSED","TRUNCATED","DSStateType","DS_STATE_",NULL};
23: const char *DSParallelTypes[] = {"REDUNDANT","SYNCHRONIZED","DISTRIBUTED","DSParallelType","DS_PARALLEL_",NULL};
24: const char *DSMatName[DS_NUM_MAT] = {"A","B","C","T","D","Q","Z","X","Y","U","V","W","E0","E1","E2","E3","E4","E5","E6","E7","E8","E9"};
25: 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};
27: /*@C
28: DSFinalizePackage - This function destroys everything in the SLEPc interface
29: to the DS package. It is called from SlepcFinalize().
31: Level: developer
33: .seealso: SlepcFinalize()
34: @*/
35: PetscErrorCode DSFinalizePackage(void)
36: {
37: PetscFunctionBegin;
38: PetscCall(PetscFunctionListDestroy(&DSList));
39: DSPackageInitialized = PETSC_FALSE;
40: DSRegisterAllCalled = PETSC_FALSE;
41: PetscFunctionReturn(PETSC_SUCCESS);
42: }
44: /*@C
45: DSInitializePackage - This function initializes everything in the DS package.
46: It is called from PetscDLLibraryRegister() when using dynamic libraries, and
47: on the first call to DSCreate() when using static libraries.
49: Level: developer
51: .seealso: SlepcInitialize()
52: @*/
53: PetscErrorCode DSInitializePackage(void)
54: {
55: char logList[256];
56: PetscBool opt,pkg;
57: PetscClassId classids[1];
59: PetscFunctionBegin;
60: if (DSPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
61: DSPackageInitialized = PETSC_TRUE;
62: /* Register Classes */
63: PetscCall(PetscClassIdRegister("Direct Solver",&DS_CLASSID));
64: /* Register Constructors */
65: PetscCall(DSRegisterAll());
66: /* Register Events */
67: PetscCall(PetscLogEventRegister("DSSolve",DS_CLASSID,&DS_Solve));
68: PetscCall(PetscLogEventRegister("DSVectors",DS_CLASSID,&DS_Vectors));
69: PetscCall(PetscLogEventRegister("DSSynchronize",DS_CLASSID,&DS_Synchronize));
70: PetscCall(PetscLogEventRegister("DSOther",DS_CLASSID,&DS_Other));
71: /* Process Info */
72: classids[0] = DS_CLASSID;
73: PetscCall(PetscInfoProcessClass("ds",1,&classids[0]));
74: /* Process summary exclusions */
75: PetscCall(PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt));
76: if (opt) {
77: PetscCall(PetscStrInList("ds",logList,',',&pkg));
78: if (pkg) PetscCall(PetscLogEventDeactivateClass(DS_CLASSID));
79: }
80: /* Register package finalizer */
81: PetscCall(PetscRegisterFinalize(DSFinalizePackage));
82: PetscFunctionReturn(PETSC_SUCCESS);
83: }
85: /*@
86: DSCreate - Creates a DS context.
88: Collective
90: Input Parameter:
91: . comm - MPI communicator
93: Output Parameter:
94: . newds - location to put the DS context
96: Level: beginner
98: Note:
99: DS objects are not intended for normal users but only for
100: advanced user that for instance implement their own solvers.
102: .seealso: DSDestroy(), DS
103: @*/
104: PetscErrorCode DSCreate(MPI_Comm comm,DS *newds)
105: {
106: DS ds;
107: PetscInt i;
109: PetscFunctionBegin;
110: PetscAssertPointer(newds,2);
111: *newds = NULL;
112: PetscCall(DSInitializePackage());
113: PetscCall(SlepcHeaderCreate(ds,DS_CLASSID,"DS","Direct Solver (or Dense System)","DS",comm,DSDestroy,DSView));
115: ds->state = DS_STATE_RAW;
116: ds->method = 0;
117: ds->compact = PETSC_FALSE;
118: ds->refined = PETSC_FALSE;
119: ds->extrarow = PETSC_FALSE;
120: ds->ld = 0;
121: ds->l = 0;
122: ds->n = 0;
123: ds->k = 0;
124: ds->t = 0;
125: ds->bs = 1;
126: ds->sc = NULL;
127: ds->pmode = DS_PARALLEL_REDUNDANT;
129: for (i=0;i<DS_NUM_MAT;i++) ds->omat[i] = NULL;
130: ds->perm = NULL;
131: ds->data = NULL;
132: ds->scset = PETSC_FALSE;
133: ds->work = NULL;
134: ds->rwork = NULL;
135: ds->iwork = NULL;
136: ds->lwork = 0;
137: ds->lrwork = 0;
138: ds->liwork = 0;
140: *newds = ds;
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
144: /*@C
145: DSSetOptionsPrefix - Sets the prefix used for searching for all
146: DS options in the database.
148: Logically Collective
150: Input Parameters:
151: + ds - the direct solver context
152: - prefix - the prefix string to prepend to all DS option requests
154: Notes:
155: A hyphen (-) must NOT be given at the beginning of the prefix name.
156: The first character of all runtime options is AUTOMATICALLY the
157: hyphen.
159: Level: advanced
161: .seealso: DSAppendOptionsPrefix()
162: @*/
163: PetscErrorCode DSSetOptionsPrefix(DS ds,const char *prefix)
164: {
165: PetscFunctionBegin;
167: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)ds,prefix));
168: PetscFunctionReturn(PETSC_SUCCESS);
169: }
171: /*@C
172: DSAppendOptionsPrefix - Appends to the prefix used for searching for all
173: DS options in the database.
175: Logically Collective
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 hyphen.
185: Level: advanced
187: .seealso: DSSetOptionsPrefix()
188: @*/
189: PetscErrorCode DSAppendOptionsPrefix(DS ds,const char *prefix)
190: {
191: PetscFunctionBegin;
193: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)ds,prefix));
194: PetscFunctionReturn(PETSC_SUCCESS);
195: }
197: /*@C
198: DSGetOptionsPrefix - Gets the prefix used for searching for all
199: DS options in the database.
201: Not Collective
203: Input Parameters:
204: . ds - the direct solver context
206: Output Parameters:
207: . prefix - pointer to the prefix string used is returned
209: Note:
210: On the Fortran side, the user should pass in a string 'prefix' of
211: sufficient length to hold the prefix.
213: Level: advanced
215: .seealso: DSSetOptionsPrefix(), DSAppendOptionsPrefix()
216: @*/
217: PetscErrorCode DSGetOptionsPrefix(DS ds,const char *prefix[])
218: {
219: PetscFunctionBegin;
221: PetscAssertPointer(prefix,2);
222: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)ds,prefix));
223: PetscFunctionReturn(PETSC_SUCCESS);
224: }
226: /*@C
227: DSSetType - Selects the type for the DS object.
229: Logically Collective
231: Input Parameters:
232: + ds - the direct solver context
233: - type - a known type
235: Level: intermediate
237: .seealso: DSGetType()
238: @*/
239: PetscErrorCode DSSetType(DS ds,DSType type)
240: {
241: PetscErrorCode (*r)(DS);
242: PetscBool match;
244: PetscFunctionBegin;
246: PetscAssertPointer(type,2);
248: PetscCall(PetscObjectTypeCompare((PetscObject)ds,type,&match));
249: if (match) PetscFunctionReturn(PETSC_SUCCESS);
251: PetscCall(PetscFunctionListFind(DSList,type,&r));
252: PetscCheck(r,PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested DS type %s",type);
254: PetscTryTypeMethod(ds,destroy);
255: PetscCall(PetscMemzero(ds->ops,sizeof(struct _DSOps)));
257: PetscCall(PetscObjectChangeTypeName((PetscObject)ds,type));
258: PetscCall((*r)(ds));
259: PetscFunctionReturn(PETSC_SUCCESS);
260: }
262: /*@C
263: DSGetType - Gets the DS type name (as a string) from the DS context.
265: Not Collective
267: Input Parameter:
268: . ds - the direct solver context
270: Output Parameter:
271: . type - name of the direct solver
273: Level: intermediate
275: .seealso: DSSetType()
276: @*/
277: PetscErrorCode DSGetType(DS ds,DSType *type)
278: {
279: PetscFunctionBegin;
281: PetscAssertPointer(type,2);
282: *type = ((PetscObject)ds)->type_name;
283: PetscFunctionReturn(PETSC_SUCCESS);
284: }
286: /*@
287: DSDuplicate - Creates a new direct solver object with the same options as
288: an existing one.
290: Collective
292: Input Parameter:
293: . ds - direct solver context
295: Output Parameter:
296: . dsnew - location to put the new DS
298: Notes:
299: DSDuplicate() DOES NOT COPY the matrices, and the new DS does not even have
300: internal arrays allocated. Use DSAllocate() to use the new DS.
302: The sorting criterion options are not copied, see DSSetSlepcSC().
304: Level: intermediate
306: .seealso: DSCreate(), DSAllocate(), DSSetSlepcSC()
307: @*/
308: PetscErrorCode DSDuplicate(DS ds,DS *dsnew)
309: {
310: PetscFunctionBegin;
312: PetscAssertPointer(dsnew,2);
313: PetscCall(DSCreate(PetscObjectComm((PetscObject)ds),dsnew));
314: if (((PetscObject)ds)->type_name) PetscCall(DSSetType(*dsnew,((PetscObject)ds)->type_name));
315: (*dsnew)->method = ds->method;
316: (*dsnew)->compact = ds->compact;
317: (*dsnew)->refined = ds->refined;
318: (*dsnew)->extrarow = ds->extrarow;
319: (*dsnew)->bs = ds->bs;
320: (*dsnew)->pmode = ds->pmode;
321: PetscFunctionReturn(PETSC_SUCCESS);
322: }
324: /*@
325: DSSetMethod - Selects the method to be used to solve the problem.
327: Logically Collective
329: Input Parameters:
330: + ds - the direct solver context
331: - meth - an index identifying the method
333: Options Database Key:
334: . -ds_method <meth> - Sets the method
336: Level: intermediate
338: .seealso: DSGetMethod()
339: @*/
340: PetscErrorCode DSSetMethod(DS ds,PetscInt meth)
341: {
342: PetscFunctionBegin;
345: PetscCheck(meth>=0,PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"The method must be a non-negative integer");
346: PetscCheck(meth<=DS_MAX_SOLVE,PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"Too large value for the method");
347: ds->method = meth;
348: PetscFunctionReturn(PETSC_SUCCESS);
349: }
351: /*@
352: DSGetMethod - Gets the method currently used in the DS.
354: Not Collective
356: Input Parameter:
357: . ds - the direct solver context
359: Output Parameter:
360: . meth - identifier of the method
362: Level: intermediate
364: .seealso: DSSetMethod()
365: @*/
366: PetscErrorCode DSGetMethod(DS ds,PetscInt *meth)
367: {
368: PetscFunctionBegin;
370: PetscAssertPointer(meth,2);
371: *meth = ds->method;
372: PetscFunctionReturn(PETSC_SUCCESS);
373: }
375: /*@
376: DSSetParallel - Selects the mode of operation in parallel runs.
378: Logically Collective
380: Input Parameters:
381: + ds - the direct solver context
382: - pmode - the parallel mode
384: Options Database Key:
385: . -ds_parallel <mode> - Sets the parallel mode, 'redundant', 'synchronized'
386: or 'distributed'
388: Notes:
389: In the 'redundant' parallel mode, all processes will make the computation
390: redundantly, starting from the same data, and producing the same result.
391: This result may be slightly different in the different processes if using a
392: multithreaded BLAS library, which may cause issues in ill-conditioned problems.
394: In the 'synchronized' parallel mode, only the first MPI process performs the
395: computation and then the computed quantities are broadcast to the other
396: processes in the communicator. This communication is not done automatically,
397: an explicit call to DSSynchronize() is required.
399: The 'distributed' parallel mode can be used in some DS types only, such
400: as the contour integral method of DSNEP. In this case, every MPI process
401: will be in charge of part of the computation.
403: Level: advanced
405: .seealso: DSSynchronize(), DSGetParallel()
406: @*/
407: PetscErrorCode DSSetParallel(DS ds,DSParallelType pmode)
408: {
409: PetscFunctionBegin;
412: ds->pmode = pmode;
413: PetscFunctionReturn(PETSC_SUCCESS);
414: }
416: /*@
417: DSGetParallel - Gets the mode of operation in parallel runs.
419: Not Collective
421: Input Parameter:
422: . ds - the direct solver context
424: Output Parameter:
425: . pmode - the parallel mode
427: Level: advanced
429: .seealso: DSSetParallel()
430: @*/
431: PetscErrorCode DSGetParallel(DS ds,DSParallelType *pmode)
432: {
433: PetscFunctionBegin;
435: PetscAssertPointer(pmode,2);
436: *pmode = ds->pmode;
437: PetscFunctionReturn(PETSC_SUCCESS);
438: }
440: /*@
441: DSSetCompact - Switch to compact storage of matrices.
443: Logically Collective
445: Input Parameters:
446: + ds - the direct solver context
447: - comp - a boolean flag
449: Notes:
450: Compact storage is used in some DS types such as DSHEP when the matrix
451: is tridiagonal. This flag can be used to indicate whether the user
452: provides the matrix entries via the compact form (the tridiagonal DS_MAT_T)
453: or the non-compact one (DS_MAT_A).
455: The default is PETSC_FALSE.
457: Level: advanced
459: .seealso: DSGetCompact()
460: @*/
461: PetscErrorCode DSSetCompact(DS ds,PetscBool comp)
462: {
463: PetscFunctionBegin;
466: ds->compact = comp;
467: PetscFunctionReturn(PETSC_SUCCESS);
468: }
470: /*@
471: DSGetCompact - Gets the compact storage flag.
473: Not Collective
475: Input Parameter:
476: . ds - the direct solver context
478: Output Parameter:
479: . comp - the flag
481: Level: advanced
483: .seealso: DSSetCompact()
484: @*/
485: PetscErrorCode DSGetCompact(DS ds,PetscBool *comp)
486: {
487: PetscFunctionBegin;
489: PetscAssertPointer(comp,2);
490: *comp = ds->compact;
491: PetscFunctionReturn(PETSC_SUCCESS);
492: }
494: /*@
495: DSSetExtraRow - Sets a flag to indicate that the matrix has one extra
496: row.
498: Logically Collective
500: Input Parameters:
501: + ds - the direct solver context
502: - ext - a boolean flag
504: Notes:
505: In Krylov methods it is useful that the matrix representing the direct solver
506: has one extra row, i.e., has dimension (n+1) x n. If this flag is activated, all
507: transformations applied to the right of the matrix also affect this additional
508: row. In that case, (n+1) must be less or equal than the leading dimension.
510: The default is PETSC_FALSE.
512: Level: advanced
514: .seealso: DSSolve(), DSAllocate(), DSGetExtraRow()
515: @*/
516: PetscErrorCode DSSetExtraRow(DS ds,PetscBool ext)
517: {
518: PetscFunctionBegin;
521: PetscCheck(!ext || ds->n==0 || ds->n!=ds->ld,PetscObjectComm((PetscObject)ds),PETSC_ERR_ORDER,"Cannot set extra row after setting n=ld");
522: ds->extrarow = ext;
523: PetscFunctionReturn(PETSC_SUCCESS);
524: }
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: {
543: PetscFunctionBegin;
545: PetscAssertPointer(ext,2);
546: *ext = ds->extrarow;
547: PetscFunctionReturn(PETSC_SUCCESS);
548: }
550: /*@
551: DSSetRefined - Sets a flag to indicate that refined vectors must be
552: computed.
554: Logically Collective
556: Input Parameters:
557: + ds - the direct solver context
558: - ref - a boolean flag
560: Notes:
561: Normally the vectors returned in DS_MAT_X are eigenvectors of the
562: projected matrix. With this flag activated, DSVectors() will return
563: the right singular vector of the smallest singular value of matrix
564: \tilde{A}-theta*I, where \tilde{A} is the extended (n+1)xn matrix
565: and theta is the Ritz value. This is used in the refined Ritz
566: approximation.
568: The default is PETSC_FALSE.
570: Level: advanced
572: .seealso: DSVectors(), DSGetRefined()
573: @*/
574: PetscErrorCode DSSetRefined(DS ds,PetscBool ref)
575: {
576: PetscFunctionBegin;
579: ds->refined = ref;
580: PetscFunctionReturn(PETSC_SUCCESS);
581: }
583: /*@
584: DSGetRefined - Gets the refined vectors flag.
586: Not Collective
588: Input Parameter:
589: . ds - the direct solver context
591: Output Parameter:
592: . ref - the flag
594: Level: advanced
596: .seealso: DSSetRefined()
597: @*/
598: PetscErrorCode DSGetRefined(DS ds,PetscBool *ref)
599: {
600: PetscFunctionBegin;
602: PetscAssertPointer(ref,2);
603: *ref = ds->refined;
604: PetscFunctionReturn(PETSC_SUCCESS);
605: }
607: /*@
608: DSSetBlockSize - Sets the block size.
610: Logically Collective
612: Input Parameters:
613: + ds - the direct solver context
614: - bs - the block size
616: Options Database Key:
617: . -ds_block_size <bs> - Sets the block size
619: Level: intermediate
621: .seealso: DSGetBlockSize()
622: @*/
623: PetscErrorCode DSSetBlockSize(DS ds,PetscInt bs)
624: {
625: PetscFunctionBegin;
628: PetscCheck(bs>0,PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"The block size must be at least one");
629: ds->bs = bs;
630: PetscFunctionReturn(PETSC_SUCCESS);
631: }
633: /*@
634: DSGetBlockSize - Gets the block size.
636: Not Collective
638: Input Parameter:
639: . ds - the direct solver context
641: Output Parameter:
642: . bs - block size
644: Level: intermediate
646: .seealso: DSSetBlockSize()
647: @*/
648: PetscErrorCode DSGetBlockSize(DS ds,PetscInt *bs)
649: {
650: PetscFunctionBegin;
652: PetscAssertPointer(bs,2);
653: *bs = ds->bs;
654: PetscFunctionReturn(PETSC_SUCCESS);
655: }
657: /*@C
658: DSSetSlepcSC - Sets the sorting criterion context.
660: Logically Collective
662: Input Parameters:
663: + ds - the direct solver context
664: - sc - a pointer to the sorting criterion context
666: Level: developer
668: .seealso: DSGetSlepcSC(), DSSort()
669: @*/
670: PetscErrorCode DSSetSlepcSC(DS ds,SlepcSC sc)
671: {
672: PetscFunctionBegin;
674: PetscAssertPointer(sc,2);
675: if (ds->sc && !ds->scset) PetscCall(PetscFree(ds->sc));
676: ds->sc = sc;
677: ds->scset = PETSC_TRUE;
678: PetscFunctionReturn(PETSC_SUCCESS);
679: }
681: /*@C
682: DSGetSlepcSC - Gets the sorting criterion context.
684: Not Collective
686: Input Parameter:
687: . ds - the direct solver context
689: Output Parameters:
690: . sc - a pointer to the sorting criterion context
692: Level: developer
694: .seealso: DSSetSlepcSC(), DSSort()
695: @*/
696: PetscErrorCode DSGetSlepcSC(DS ds,SlepcSC *sc)
697: {
698: PetscFunctionBegin;
700: PetscAssertPointer(sc,2);
701: if (!ds->sc) PetscCall(PetscNew(&ds->sc));
702: *sc = ds->sc;
703: PetscFunctionReturn(PETSC_SUCCESS);
704: }
706: /*@
707: DSSetFromOptions - Sets DS options from the options database.
709: Collective
711: Input Parameters:
712: . ds - the direct solver context
714: Notes:
715: To see all options, run your program with the -help option.
717: Level: beginner
719: .seealso: DSSetOptionsPrefix()
720: @*/
721: PetscErrorCode DSSetFromOptions(DS ds)
722: {
723: PetscInt bs,meth;
724: PetscBool flag;
725: DSParallelType pmode;
727: PetscFunctionBegin;
729: PetscCall(DSRegisterAll());
730: /* Set default type (we do not allow changing it with -ds_type) */
731: if (!((PetscObject)ds)->type_name) PetscCall(DSSetType(ds,DSNHEP));
732: PetscObjectOptionsBegin((PetscObject)ds);
734: PetscCall(PetscOptionsInt("-ds_block_size","Block size for the dense system solver","DSSetBlockSize",ds->bs,&bs,&flag));
735: if (flag) PetscCall(DSSetBlockSize(ds,bs));
737: PetscCall(PetscOptionsInt("-ds_method","Method to be used for the dense system","DSSetMethod",ds->method,&meth,&flag));
738: if (flag) PetscCall(DSSetMethod(ds,meth));
740: PetscCall(PetscOptionsEnum("-ds_parallel","Operation mode in parallel runs","DSSetParallel",DSParallelTypes,(PetscEnum)ds->pmode,(PetscEnum*)&pmode,&flag));
741: if (flag) PetscCall(DSSetParallel(ds,pmode));
743: PetscTryTypeMethod(ds,setfromoptions,PetscOptionsObject);
744: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)ds,PetscOptionsObject));
745: PetscOptionsEnd();
746: PetscFunctionReturn(PETSC_SUCCESS);
747: }
749: /*@C
750: DSView - Prints the DS data structure.
752: Collective
754: Input Parameters:
755: + ds - the direct solver context
756: - viewer - optional visualization context
758: Note:
759: The available visualization contexts include
760: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
761: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
762: output where only the first processor opens
763: the file. All other processors send their
764: data to the first processor to print.
766: The user can open an alternative visualization context with
767: PetscViewerASCIIOpen() - output to a specified file.
769: Level: beginner
771: .seealso: DSViewMat()
772: @*/
773: PetscErrorCode DSView(DS ds,PetscViewer viewer)
774: {
775: PetscBool isascii;
776: PetscViewerFormat format;
777: PetscMPIInt size;
779: PetscFunctionBegin;
781: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)ds),&viewer));
783: PetscCheckSameComm(ds,1,viewer,2);
784: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
785: if (isascii) {
786: PetscCall(PetscViewerGetFormat(viewer,&format));
787: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)ds,viewer));
788: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ds),&size));
789: if (size>1) PetscCall(PetscViewerASCIIPrintf(viewer," parallel operation mode: %s\n",DSParallelTypes[ds->pmode]));
790: if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
791: PetscCall(PetscViewerASCIIPrintf(viewer," current state: %s\n",DSStateTypes[ds->state]));
792: PetscCall(PetscViewerASCIIPrintf(viewer," dimensions: ld=%" PetscInt_FMT ", n=%" PetscInt_FMT ", l=%" PetscInt_FMT ", k=%" PetscInt_FMT,ds->ld,ds->n,ds->l,ds->k));
793: if (ds->state==DS_STATE_TRUNCATED) PetscCall(PetscViewerASCIIPrintf(viewer,", t=%" PetscInt_FMT "\n",ds->t));
794: else PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
795: PetscCall(PetscViewerASCIIPrintf(viewer," flags:%s%s%s\n",ds->compact?" compact":"",ds->extrarow?" extrarow":"",ds->refined?" refined":""));
796: }
797: PetscCall(PetscViewerASCIIPushTab(viewer));
798: PetscTryTypeMethod(ds,view,viewer);
799: PetscCall(PetscViewerASCIIPopTab(viewer));
800: }
801: PetscFunctionReturn(PETSC_SUCCESS);
802: }
804: /*@C
805: DSViewFromOptions - View from options
807: Collective
809: Input Parameters:
810: + ds - the direct solver context
811: . obj - optional object
812: - name - command line option
814: Level: intermediate
816: .seealso: DSView(), DSCreate()
817: @*/
818: PetscErrorCode DSViewFromOptions(DS ds,PetscObject obj,const char name[])
819: {
820: PetscFunctionBegin;
822: PetscCall(PetscObjectViewFromOptions((PetscObject)ds,obj,name));
823: PetscFunctionReturn(PETSC_SUCCESS);
824: }
826: /*@
827: DSAllocate - Allocates memory for internal storage or matrices in DS.
829: Logically Collective
831: Input Parameters:
832: + ds - the direct solver context
833: - ld - leading dimension (maximum allowed dimension for the matrices, including
834: the extra row if present)
836: Note:
837: If the leading dimension is different from a previously set value, then
838: all matrices are destroyed with DSReset().
840: Level: intermediate
842: .seealso: DSGetLeadingDimension(), DSSetDimensions(), DSSetExtraRow(), DSReset()
843: @*/
844: PetscErrorCode DSAllocate(DS ds,PetscInt ld)
845: {
846: PetscFunctionBegin;
850: PetscCheck(ld>0,PetscObjectComm((PetscObject)ds),PETSC_ERR_ARG_OUTOFRANGE,"Leading dimension should be at least one");
851: if (ld!=ds->ld) {
852: PetscCall(PetscInfo(ds,"Allocating memory with leading dimension=%" PetscInt_FMT "\n",ld));
853: PetscCall(DSReset(ds));
854: ds->ld = ld;
855: PetscUseTypeMethod(ds,allocate,ld);
856: }
857: PetscFunctionReturn(PETSC_SUCCESS);
858: }
860: /*@
861: DSReset - Resets the DS context to the initial state.
863: Collective
865: Input Parameter:
866: . ds - the direct solver context
868: Note:
869: All data structures with size depending on the leading dimension
870: of DSAllocate() are released.
872: Level: advanced
874: .seealso: DSDestroy(), DSAllocate()
875: @*/
876: PetscErrorCode DSReset(DS ds)
877: {
878: PetscInt i;
880: PetscFunctionBegin;
882: if (!ds) PetscFunctionReturn(PETSC_SUCCESS);
883: ds->state = DS_STATE_RAW;
884: ds->ld = 0;
885: ds->l = 0;
886: ds->n = 0;
887: ds->k = 0;
888: for (i=0;i<DS_NUM_MAT;i++) PetscCall(MatDestroy(&ds->omat[i]));
889: PetscCall(PetscFree(ds->perm));
890: PetscFunctionReturn(PETSC_SUCCESS);
891: }
893: /*@C
894: DSDestroy - Destroys DS context that was created with DSCreate().
896: Collective
898: Input Parameter:
899: . ds - the direct solver context
901: Level: beginner
903: .seealso: DSCreate()
904: @*/
905: PetscErrorCode DSDestroy(DS *ds)
906: {
907: PetscFunctionBegin;
908: if (!*ds) PetscFunctionReturn(PETSC_SUCCESS);
910: if (--((PetscObject)*ds)->refct > 0) { *ds = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
911: PetscCall(DSReset(*ds));
912: PetscTryTypeMethod(*ds,destroy);
913: PetscCall(PetscFree((*ds)->work));
914: PetscCall(PetscFree((*ds)->rwork));
915: PetscCall(PetscFree((*ds)->iwork));
916: if (!(*ds)->scset) PetscCall(PetscFree((*ds)->sc));
917: PetscCall(PetscHeaderDestroy(ds));
918: PetscFunctionReturn(PETSC_SUCCESS);
919: }
921: /*@C
922: DSRegister - Adds a direct solver to the DS package.
924: Not Collective
926: Input Parameters:
927: + name - name of a new user-defined DS
928: - function - routine to create context
930: Notes:
931: DSRegister() may be called multiple times to add several user-defined
932: direct solvers.
934: Level: advanced
936: .seealso: DSRegisterAll()
937: @*/
938: PetscErrorCode DSRegister(const char *name,PetscErrorCode (*function)(DS))
939: {
940: PetscFunctionBegin;
941: PetscCall(DSInitializePackage());
942: PetscCall(PetscFunctionListAdd(&DSList,name,function));
943: PetscFunctionReturn(PETSC_SUCCESS);
944: }
946: SLEPC_EXTERN PetscErrorCode DSCreate_HEP(DS);
947: SLEPC_EXTERN PetscErrorCode DSCreate_NHEP(DS);
948: SLEPC_EXTERN PetscErrorCode DSCreate_GHEP(DS);
949: SLEPC_EXTERN PetscErrorCode DSCreate_GHIEP(DS);
950: SLEPC_EXTERN PetscErrorCode DSCreate_GNHEP(DS);
951: SLEPC_EXTERN PetscErrorCode DSCreate_NHEPTS(DS);
952: SLEPC_EXTERN PetscErrorCode DSCreate_SVD(DS);
953: SLEPC_EXTERN PetscErrorCode DSCreate_HSVD(DS);
954: SLEPC_EXTERN PetscErrorCode DSCreate_GSVD(DS);
955: SLEPC_EXTERN PetscErrorCode DSCreate_PEP(DS);
956: SLEPC_EXTERN PetscErrorCode DSCreate_NEP(DS);
958: /*@C
959: DSRegisterAll - Registers all of the direct solvers in the DS package.
961: Not Collective
963: Level: advanced
965: .seealso: DSRegister()
966: @*/
967: PetscErrorCode DSRegisterAll(void)
968: {
969: PetscFunctionBegin;
970: if (DSRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS);
971: DSRegisterAllCalled = PETSC_TRUE;
972: PetscCall(DSRegister(DSHEP,DSCreate_HEP));
973: PetscCall(DSRegister(DSNHEP,DSCreate_NHEP));
974: PetscCall(DSRegister(DSGHEP,DSCreate_GHEP));
975: PetscCall(DSRegister(DSGHIEP,DSCreate_GHIEP));
976: PetscCall(DSRegister(DSGNHEP,DSCreate_GNHEP));
977: PetscCall(DSRegister(DSNHEPTS,DSCreate_NHEPTS));
978: PetscCall(DSRegister(DSSVD,DSCreate_SVD));
979: PetscCall(DSRegister(DSHSVD,DSCreate_HSVD));
980: PetscCall(DSRegister(DSGSVD,DSCreate_GSVD));
981: PetscCall(DSRegister(DSPEP,DSCreate_PEP));
982: PetscCall(DSRegister(DSNEP,DSCreate_NEP));
983: PetscFunctionReturn(PETSC_SUCCESS);
984: }