Actual source code: sles.c
1: /*$Id: sles.c,v 1.151 2001/08/07 03:03:26 balay Exp $*/
3: #include src/sles/slesimpl.h
5: /* Logging support */
6: int SLES_COOKIE;
7: int SLES_SetUp, SLES_Solve;
9: static int SLESPublish_Petsc(PetscObject obj)
10: {
11: #if defined(PETSC_HAVE_AMS)
12: int ierr;
13: #endif
14:
16: #if defined(PETSC_HAVE_AMS)
17: PetscObjectPublishBaseBegin(obj);
18: PetscObjectPublishBaseEnd(obj);
19: #endif
20: return(0);
21: }
23: /*@C
24: SLESView - Prints the SLES data structure.
26: Collective on SLES
28: Input Parameters:
29: + SLES - the SLES context
30: - viewer - optional visualization context
32: Options Database Key:
33: . -sles_view - Calls SLESView() at end of SLESSolve()
35: Note:
36: The available visualization contexts include
37: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
38: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
39: output where only the first processor opens
40: the file. All other processors send their
41: data to the first processor to print.
43: The user can open alternative visualization contexts with
44: . PetscViewerASCIIOpen() - output to a specified file
46: Level: beginner
48: .keywords: SLES, view
50: .seealso: PetscViewerASCIIOpen()
51: @*/
52: int SLESView(SLES sles,PetscViewer viewer)
53: {
54: KSP ksp;
55: PC pc;
56: int ierr;
60: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(sles->comm);
63: SLESGetPC(sles,&pc);
64: SLESGetKSP(sles,&ksp);
65: KSPView(ksp,viewer);
66: PCView(pc,viewer);
67: return(0);
68: }
70: /*@C
71: SLESSetOptionsPrefix - Sets the prefix used for searching for all
72: SLES options in the database.
74: Collective on SLES
76: Input Parameter:
77: + sles - the SLES context
78: - prefix - the prefix to prepend to all option names
80: Notes:
81: A hyphen (-) must NOT be given at the beginning of the prefix name.
82: The first character of all runtime options is AUTOMATICALLY the
83: hyphen.
85: This prefix is particularly useful for nested use of SLES. For
86: example, the block Jacobi and block diagonal preconditioners use
87: the prefix "sub_" for options relating to the individual blocks.
89: Level: intermediate
91: .keywords: SLES, set, options, prefix, database
93: .seealso: SLESAppendOptionsPrefix(), SLESGetOptionsPrefix()
94: @*/
95: int SLESSetOptionsPrefix(SLES sles,char *prefix)
96: {
101: PetscObjectSetOptionsPrefix((PetscObject)sles,prefix);
102: KSPSetOptionsPrefix(sles->ksp,prefix);
103: PCSetOptionsPrefix(sles->pc,prefix);
104: return(0);
105: }
107: /*@C
108: SLESAppendOptionsPrefix - Appends to the prefix used for searching for all
109: SLES options in the database.
111: Collective on SLES
113: Input Parameter:
114: + sles - the SLES context
115: - prefix - the prefix to prepend to all option names
117: Notes:
118: A hyphen (-) must NOT be given at the beginning of the prefix name.
119: The first character of all runtime options is AUTOMATICALLY the
120: hyphen.
122: This prefix is particularly useful for nested use of SLES. For
123: example, the block Jacobi and block diagonal preconditioners use
124: the prefix "sub_" for options relating to the individual blocks.
126: Level: intermediate
128: .keywords: SLES, append, options, prefix, database
130: .seealso: SLESSetOptionsPrefix(), SLESGetOptionsPrefix()
131: @*/
132: int SLESAppendOptionsPrefix(SLES sles,char *prefix)
133: {
138: PetscObjectAppendOptionsPrefix((PetscObject)sles,prefix);
139: KSPAppendOptionsPrefix(sles->ksp,prefix);
140: PCAppendOptionsPrefix(sles->pc,prefix);
141: return(0);
142: }
144: /*@C
145: SLESGetOptionsPrefix - Gets the prefix used for searching for all
146: SLES options in the database.
148: Not Collective
150: Input Parameter:
151: . sles - the SLES context
153: Output Parameter:
154: . prefix - pointer to the prefix string used
156: Notes:
157: This prefix is particularly useful for nested use of SLES. For
158: example, the block Jacobi and block diagonal preconditioners use
159: the prefix "sub" for options relating to the individual blocks.
161: On the fortran side, the user should pass in a string 'prifix' of
162: sufficient length to hold the prefix.
164: Level: intermediate
166: .keywords: SLES, get, options, prefix, database
168: .seealso: SLESSetOptionsPrefix()
169: @*/
170: int SLESGetOptionsPrefix(SLES sles,char **prefix)
171: {
176: PetscObjectGetOptionsPrefix((PetscObject)sles,prefix);
177: return(0);
178: }
180: /*@
181: SLESSetFromOptions - Sets various SLES parameters from user options.
182: Also takes all KSP and PC options.
184: Collective on SLES
186: Input Parameter:
187: . sles - the SLES context
189: Level: beginner
191: .keywords: SLES, set, options, database
193: .seealso: KSPSetFromOptions(),
194: PCSetFromOptions(), SLESGetPC(), SLESGetKSP()
195: @*/
196: int SLESSetFromOptions(SLES sles)
197: {
198: int ierr;
199: PetscTruth flag;
203: PetscOptionsBegin(sles->comm,sles->prefix,"Linear solver (SLES) options","SLES");
204: PetscOptionsName("-sles_diagonal_scale","Diagonal scale matrix before building preconditioner","SLESSetDiagonalScale",&flag);
205: if (flag) {
206: SLESSetDiagonalScale(sles,PETSC_TRUE);
207: }
208: PetscOptionsName("-sles_diagonal_scale_fix","Fix diagonaled scaled matrix after solve","SLESSetDiagonalScaleFix",&flag);
209: if (flag) {
210: SLESSetDiagonalScaleFix(sles);
211: }
212: PetscOptionsEnd();
213: KSPSetPC(sles->ksp,sles->pc);
214: KSPSetFromOptions(sles->ksp);
215: PCSetFromOptions(sles->pc);
216: return(0);
217: }
219: /*@C
220: SLESCreate - Creates a linear equation solver context.
222: Collective on MPI_Comm
224: Input Parameter:
225: . comm - MPI communicator
227: Output Parameter:
228: . sles - the newly created SLES context
230: Level: beginner
232: .keywords: SLES, create, context
234: .seealso: SLESSolve(), SLESDestroy(), SLES
235: @*/
236: int SLESCreate(MPI_Comm comm,SLES *outsles)
237: {
238: int ierr;
239: SLES sles;
243: *outsles = 0;
244: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
245: SLESInitializePackage(PETSC_NULL);
246: #endif
248: PetscHeaderCreate(sles,_p_SLES,int,SLES_COOKIE,0,"SLES",comm,SLESDestroy,SLESView);
249: PetscLogObjectCreate(sles);
250: sles->bops->publish = SLESPublish_Petsc;
251: KSPCreate(comm,&sles->ksp);
252: PCCreate(comm,&sles->pc);
253: PetscLogObjectParent(sles,sles->ksp);
254: PetscLogObjectParent(sles,sles->pc);
255: sles->setupcalled = 0;
256: sles->dscale = PETSC_FALSE;
257: sles->dscalefix = PETSC_FALSE;
258: sles->diagonal = PETSC_NULL;
259: *outsles = sles;
260: PetscPublishAll(sles);
261: return(0);
262: }
264: /*@C
265: SLESDestroy - Destroys the SLES context.
267: Collective on SLES
269: Input Parameters:
270: . sles - the SLES context
272: Level: beginner
274: .keywords: SLES, destroy, context
276: .seealso: SLESCreate(), SLESSolve()
277: @*/
278: int SLESDestroy(SLES sles)
279: {
284: if (--sles->refct > 0) return(0);
286: PetscObjectDepublish(sles);
287: KSPDestroy(sles->ksp);
288: PCDestroy(sles->pc);
289: if (sles->diagonal) {VecDestroy(sles->diagonal);}
290: PetscLogObjectDestroy(sles);
291: PetscHeaderDestroy(sles);
292: return(0);
293: }
295: /*@
296: SLESSetUp - Performs set up required for solving a linear system.
298: Collective on SLES
300: Input Parameters:
301: + sles - the SLES context
302: . b - the right hand side
303: - x - location to hold solution
305: Note:
306: For basic use of the SLES solvers the user need not explicitly call
307: SLESSetUp(), since these actions will automatically occur during
308: the call to SLESSolve(). However, if one wishes to generate
309: performance data for this computational phase (for example, for
310: incomplete factorization using the ILU preconditioner) using the
311: PETSc log facilities, calling SLESSetUp() is required.
313: Level: advanced
315: .keywords: SLES, solve, linear system
317: .seealso: SLESCreate(), SLESDestroy(), SLESDestroy()
318: @*/
319: int SLESSetUp(SLES sles,Vec b,Vec x)
320: {
322: KSP ksp;
323: PC pc;
324: Mat mat,pmat;
333: ksp = sles->ksp; pc = sles->pc;
334: KSPSetRhs(ksp,b);
335: KSPSetSolution(ksp,x);
336: if (!sles->setupcalled) {
337: PetscLogEventBegin(SLES_SetUp,sles,b,x,0);
338: KSPSetPC(ksp,pc);
339: PCSetVector(pc,b);
340: KSPSetUp(sles->ksp);
342: /* scale the matrix if requested */
343: if (sles->dscale) {
344: PCGetOperators(pc,&mat,&pmat,PETSC_NULL);
345: if (mat == pmat) {
346: PetscScalar *xx;
347: int i,n;
348: PetscTruth zeroflag = PETSC_FALSE;
351: if (!sles->diagonal) { /* allocate vector to hold diagonal */
352: VecDuplicate(b,&sles->diagonal);
353: }
354: MatGetDiagonal(mat,sles->diagonal);
355: VecGetLocalSize(sles->diagonal,&n);
356: VecGetArray(sles->diagonal,&xx);
357: for (i=0; i<n; i++) {
358: if (xx[i] != 0.0) xx[i] = 1.0/sqrt(PetscAbsScalar(xx[i]));
359: else {
360: xx[i] = 1.0;
361: zeroflag = PETSC_TRUE;
362: }
363: }
364: VecRestoreArray(sles->diagonal,&xx);
365: if (zeroflag) {
366: PetscLogInfo(pc,"SLESSetUp:Zero detected in diagonal of matrix, using 1 at those locationsn");
367: }
368: MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
369: sles->dscalefix2 = PETSC_FALSE;
370: }
371: }
373: PCSetUp(sles->pc);
374: sles->setupcalled = 1;
375: PetscLogEventEnd(SLES_SetUp,sles,b,x,0);
376: }
377: return(0);
378: }
380: /*@
381: SLESSolve - Solves a linear system.
383: Collective on SLES
385: Input Parameters:
386: + sles - the SLES context
387: - b - the right hand side
389: Output Parameters:
390: + x - the approximate solution
391: - its - the number of iterations until termination
393: Options Database:
394: + -sles_diagonal_scale - diagonally scale linear system before solving
395: . -sles_diagonal_scale_fix - unscale the matrix and right hand side when done
396: . -sles_view - print information about the solver used
397: - -sles_view_binary - save the matrix and right hand side to the default binary file (can be
398: read later with src/sles/examples/tutorials/ex10.c for testing solvers)
400: Notes:
401: Call KSPGetConvergedReason() to determine if the solver converged or failed and
402: why.
404: On return, the parameter "its" contains the iteration number at which convergence
405: was successfully reached or divergence/failure was detected
407: If using a direct method (e.g., via the KSP solver
408: KSPPREONLY and a preconditioner such as PCLU/PCILU),
409: then its=1. See KSPSetTolerances() and KSPDefaultConverged()
410: for more details. Also, see KSPSolve() for more details
411: about iterative solver options.
413: Setting a Nonzero Initial Guess:
414: By default, SLES assumes an initial guess of zero by zeroing
415: the initial value for the solution vector, x. To use a nonzero
416: initial guess, the user must call
417: .vb
418: SLESGetKSP(sles,&ksp);
419: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
420: .ve
422: Solving Successive Linear Systems:
423: When solving multiple linear systems of the same size with the
424: same method, several options are available.
426: (1) To solve successive linear systems having the SAME
427: preconditioner matrix (i.e., the same data structure with exactly
428: the same matrix elements) but different right-hand-side vectors,
429: the user should simply call SLESSolve() multiple times. The
430: preconditioner setup operations (e.g., factorization for ILU) will
431: be done during the first call to SLESSolve() only; such operations
432: will NOT be repeated for successive solves.
434: (2) To solve successive linear systems that have DIFFERENT
435: preconditioner matrices (i.e., the matrix elements and/or the
436: matrix data structure change), the user MUST call SLESSetOperators()
437: and SLESSolve() for each solve. See SLESSetOperators() for
438: options that can save work for such cases.
440: Level: beginner
442: .keywords: SLES, solve, linear system
444: .seealso: SLESCreate(), SLESSetOperators(), SLESGetKSP(), KSPSetTolerances(),
445: KSPDefaultConverged(), KSPSetInitialGuessNonzero(), KSPGetResidualNorm(),
446: KSPSolve()
447: @*/
448: int SLESSolve(SLES sles,Vec b,Vec x,int *its)
449: {
450: int ierr;
451: PetscTruth flg;
452: KSP ksp;
453: PC pc;
461: if (b == x) SETERRQ(PETSC_ERR_ARG_IDN,"b and x must be different vectors");
463: PetscOptionsHasName(sles->prefix,"-sles_view_binary",&flg);
464: if (flg) {
465: Mat mat;
466: PCGetOperators(sles->pc,&mat,PETSC_NULL,PETSC_NULL);
467: MatView(mat,PETSC_VIEWER_BINARY_(sles->comm));
468: VecView(b,PETSC_VIEWER_BINARY_(sles->comm));
469: }
470: ksp = sles->ksp;
471: pc = sles->pc;
472: SLESSetUp(sles,b,x);
473: SLESSetUpOnBlocks(sles);
474: PetscLogEventBegin(SLES_Solve,sles,b,x,0);
475: /* diagonal scale RHS if called for */
476: if (sles->dscale) {
477: VecPointwiseMult(sles->diagonal,b,b);
478: /* second time in, but matrix was scaled back to original */
479: if (sles->dscalefix && sles->dscalefix2) {
480: Mat mat;
482: PCGetOperators(pc,&mat,PETSC_NULL,PETSC_NULL);
483: MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
484: }
485: }
486: PCPreSolve(pc,ksp);
487: KSPSolve(ksp,its);
488: PCPostSolve(pc,ksp);
489: /* diagonal scale solution if called for */
490: if (sles->dscale) {
491: VecPointwiseMult(sles->diagonal,x,x);
492: /* unscale right hand side and matrix */
493: if (sles->dscalefix) {
494: Mat mat;
496: VecReciprocal(sles->diagonal);
497: VecPointwiseMult(sles->diagonal,b,b);
498: PCGetOperators(pc,&mat,PETSC_NULL,PETSC_NULL);
499: MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
500: VecReciprocal(sles->diagonal);
501: sles->dscalefix2 = PETSC_TRUE;
502: }
503: }
504: PetscLogEventEnd(SLES_Solve,sles,b,x,0);
505: PetscOptionsHasName(sles->prefix,"-sles_view",&flg);
506: if (flg && !PetscPreLoadingOn) { SLESView(sles,PETSC_VIEWER_STDOUT_WORLD); }
507: return(0);
508: }
510: /*@
511: SLESSolveTranspose - Solves the transpose of a linear system.
513: Collective on SLES
515: Input Parameters:
516: + sles - the SLES context
517: - b - the right hand side
519: Output Parameters:
520: + x - the approximate solution
521: - its - the number of iterations until termination
523: Notes:
524: On return, the parameter "its" contains
525: + <its> - the iteration number at which convergence
526: was successfully reached, or
527: - <-its> - the negative of the iteration at which
528: divergence or breakdown was detected.
530: Currently works only with the KSPPREONLY method.
532: Setting a Nonzero Initial Guess:
533: By default, SLES assumes an initial guess of zero by zeroing
534: the initial value for the solution vector, x. To use a nonzero
535: initial guess, the user must call
536: .vb
537: SLESGetKSP(sles,&ksp);
538: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
539: .ve
541: Solving Successive Linear Systems:
542: When solving multiple linear systems of the same size with the
543: same method, several options are available.
545: (1) To solve successive linear systems having the SAME
546: preconditioner matrix (i.e., the same data structure with exactly
547: the same matrix elements) but different right-hand-side vectors,
548: the user should simply call SLESSolve() multiple times. The
549: preconditioner setup operations (e.g., factorization for ILU) will
550: be done during the first call to SLESSolve() only; such operations
551: will NOT be repeated for successive solves.
553: (2) To solve successive linear systems that have DIFFERENT
554: preconditioner matrices (i.e., the matrix elements and/or the
555: matrix data structure change), the user MUST call SLESSetOperators()
556: and SLESSolve() for each solve. See SLESSetOperators() for
557: options that can save work for such cases.
559: Level: intermediate
561: .keywords: SLES, solve, linear system
563: .seealso: SLESCreate(), SLESSetOperators(), SLESGetKSP(), KSPSetTolerances(),
564: KSPDefaultConverged(), KSPSetInitialGuessNonzero(), KSPGetResidualNorm(),
565: KSPSolve()
566: @*/
567: int SLESSolveTranspose(SLES sles,Vec b,Vec x,int *its)
568: {
569: int ierr;
570: PetscTruth flg;
571: KSP ksp;
572: PC pc;
578: if (b == x) SETERRQ(PETSC_ERR_ARG_IDN,"b and x must be different vectors");
582: ksp = sles->ksp; pc = sles->pc;
583: KSPSetRhs(ksp,b);
584: KSPSetSolution(ksp,x);
585: KSPSetPC(ksp,pc);
586: if (!sles->setupcalled) {
587: SLESSetUp(sles,b,x);
588: }
589: PCSetUpOnBlocks(pc);
590: PetscLogEventBegin(SLES_Solve,sles,b,x,0);
591: KSPSolveTranspose(ksp,its);
592: PetscLogEventEnd(SLES_Solve,sles,b,x,0);
593: PetscOptionsHasName(sles->prefix,"-sles_view",&flg);
594: if (flg) { SLESView(sles,PETSC_VIEWER_STDOUT_WORLD); }
595: return(0);
596: }
598: /*@C
599: SLESGetKSP - Returns the KSP context for a SLES solver.
601: Not Collective, but if KSP will be a parallel object if SLES is
603: Input Parameter:
604: . sles - the SLES context
606: Output Parameter:
607: . ksp - the Krylov space context
609: Notes:
610: The user can then directly manipulate the KSP context to set various
611: options (e.g., by calling KSPSetType()), etc.
613: Level: beginner
614:
615: .keywords: SLES, get, KSP, context
617: .seealso: SLESGetPC()
618: @*/
619: int SLESGetKSP(SLES sles,KSP *ksp)
620: {
623: *ksp = sles->ksp;
624: return(0);
625: }
627: /*@C
628: SLESGetPC - Returns the preconditioner (PC) context for a SLES solver.
630: Not Collective, but if PC will be a parallel object if SLES is
632: Input Parameter:
633: . sles - the SLES context
635: Output Parameter:
636: . pc - the preconditioner context
638: Notes:
639: The user can then directly manipulate the PC context to set various
640: options (e.g., by calling PCSetType()), etc.
642: Level: beginner
644: .keywords: SLES, get, PC, context
646: .seealso: SLESGetKSP()
647: @*/
648: int SLESGetPC(SLES sles,PC *pc)
649: {
652: *pc = sles->pc;
653: return(0);
654: }
656: #include src/mat/matimpl.h
657: /*@
658: SLESSetOperators - Sets the matrix associated with the linear system
659: and a (possibly) different one associated with the preconditioner.
661: Collective on SLES and Mat
663: Input Parameters:
664: + sles - the SLES context
665: . Amat - the matrix associated with the linear system
666: . Pmat - the matrix to be used in constructing the preconditioner, usually the
667: same as Amat.
668: - flag - flag indicating information about the preconditioner matrix structure
669: during successive linear solves. This flag is ignored the first time a
670: linear system is solved, and thus is irrelevant when solving just one linear
671: system.
673: Notes:
674: The flag can be used to eliminate unnecessary work in the preconditioner
675: during the repeated solution of linear systems of the same size. The
676: available options are
677: $ SAME_PRECONDITIONER -
678: $ Pmat is identical during successive linear solves.
679: $ This option is intended for folks who are using
680: $ different Amat and Pmat matrices and want to reuse the
681: $ same preconditioner matrix. For example, this option
682: $ saves work by not recomputing incomplete factorization
683: $ for ILU/ICC preconditioners.
684: $ SAME_NONZERO_PATTERN -
685: $ Pmat has the same nonzero structure during
686: $ successive linear solves.
687: $ DIFFERENT_NONZERO_PATTERN -
688: $ Pmat does not have the same nonzero structure.
690: Caution:
691: If you specify SAME_NONZERO_PATTERN, PETSc believes your assertion
692: and does not check the structure of the matrix. If you erroneously
693: claim that the structure is the same when it actually is not, the new
694: preconditioner will not function correctly. Thus, use this optimization
695: feature carefully!
697: If in doubt about whether your preconditioner matrix has changed
698: structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
700: Level: beginner
702: .keywords: SLES, set, operators, matrix, preconditioner, linear system
704: .seealso: SLESSolve(), SLESGetPC(), PCGetOperators()
705: @*/
706: int SLESSetOperators(SLES sles,Mat Amat,Mat Pmat,MatStructure flag)
707: {
712: PCSetOperators(sles->pc,Amat,Pmat,flag);
713: sles->setupcalled = 0; /* so that next solve call will call setup */
714: return(0);
715: }
717: /*@
718: SLESSetUpOnBlocks - Sets up the preconditioner for each block in
719: the block Jacobi, block Gauss-Seidel, and overlapping Schwarz
720: methods.
722: Collective on SLES
724: Input Parameter:
725: . sles - the SLES context
727: Notes:
728: SLESSetUpOnBlocks() is a routine that the user can optinally call for
729: more precise profiling (via -log_summary) of the setup phase for these
730: block preconditioners. If the user does not call SLESSetUpOnBlocks(),
731: it will automatically be called from within SLESSolve().
732:
733: Calling SLESSetUpOnBlocks() is the same as calling PCSetUpOnBlocks()
734: on the PC context within the SLES context.
736: Level: advanced
738: .keywords: SLES, setup, blocks
740: .seealso: PCSetUpOnBlocks(), SLESSetUp(), PCSetUp()
741: @*/
742: int SLESSetUpOnBlocks(SLES sles)
743: {
745: PC pc;
749: SLESGetPC(sles,&pc);
750: PCSetUpOnBlocks(pc);
751: return(0);
752: }
754: /*@C
755: SLESSetDiagonalScale - Tells SLES to diagonally scale the system
756: before solving. This actually CHANGES the matrix (and right hand side).
758: Collective on SLES
760: Input Parameter:
761: + sles - the SLES context
762: - scale - PETSC_TRUE or PETSC_FALSE
764: Notes:
765: BE CAREFUL with this routine: it actually scales the matrix and right
766: hand side that define the system. After the system is solved the matrix
767: and right hand side remain scaled.
769: This routine is only used if the matrix and preconditioner matrix are
770: the same thing.
771:
772: If you use this with the PCType Eisenstat preconditioner than you can
773: use the PCEisenstatNoDiagonalScaling() option, or -pc_eisenstat_no_diagonal_scaling
774: to save some unneeded, redundant flops.
776: Level: intermediate
778: .keywords: SLES, set, options, prefix, database
780: .seealso: SLESGetDiagonalScale(), SLESSetDiagonalScaleFix()
781: @*/
782: int SLESSetDiagonalScale(SLES sles,PetscTruth scale)
783: {
786: sles->dscale = scale;
787: return(0);
788: }
790: /*@C
791: SLESGetDiagonalScale - Checks if SLES solver scales the matrix and
792: right hand side
794: Not Collective
796: Input Parameter:
797: . sles - the SLES context
799: Output Parameter:
800: . scale - PETSC_TRUE or PETSC_FALSE
802: Notes:
803: BE CAREFUL with this routine: it actually scales the matrix and right
804: hand side that define the system. After the system is solved the matrix
805: and right hand side remain scaled.
807: This routine is only used if the matrix and preconditioner matrix are
808: the same thing.
810: Level: intermediate
812: .keywords: SLES, set, options, prefix, database
814: .seealso: SLESSetDiagonalScale(), SLESSetDiagonalScaleFix()
815: @*/
816: int SLESGetDiagonalScale(SLES sles,PetscTruth *scale)
817: {
820: *scale = sles->dscale;
821: return(0);
822: }
824: /*@C
825: SLESSetDiagonalScaleFix - Tells SLES to diagonally scale the system
826: back after solving.
828: Collective on SLES
830: Input Parameter:
831: . sles - the SLES context
834: Notes:
835: Must be called after SLESSetDiagonalScale()
837: Using this will slow things down, because it rescales the matrix before and
838: after each linear solve. This is intended mainly for testing to allow one
839: to easily get back the original system to make sure the solution computed is
840: accurate enough.
842: This routine is only used if the matrix and preconditioner matrix are
843: the same thing.
845: Level: intermediate
847: .keywords: SLES, set, options, prefix, database
849: .seealso: SLESGetDiagonalScale(), SLESSetDiagonalScale()
850: @*/
851: int SLESSetDiagonalScaleFix(SLES sles)
852: {
855: if (!sles->dscale) {
856: SETERRQ(1,"Must call after SLESSetDiagonalScale(); if setting with command line options you must setn
857: -sles_diagonal_scale option if you set the -sles_diagonal_scale_fix option");
858: }
859: sles->dscalefix = PETSC_TRUE;
860: return(0);
861: }