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;
11: static int SLESPublish_Petsc(PetscObject obj)
12: {
13: #if defined(PETSC_HAVE_AMS)
14: int ierr;
15: #endif
16:
18: #if defined(PETSC_HAVE_AMS)
19: PetscObjectPublishBaseBegin(obj);
20: PetscObjectPublishBaseEnd(obj);
21: #endif
22: return(0);
23: }
27: /*@C
28: SLESView - Prints the SLES data structure.
30: Collective on SLES
32: Input Parameters:
33: + SLES - the SLES context
34: - viewer - optional visualization context
36: Options Database Key:
37: . -sles_view - Calls SLESView() at end of SLESSolve()
39: Note:
40: The available visualization contexts include
41: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
42: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
43: output where only the first processor opens
44: the file. All other processors send their
45: data to the first processor to print.
47: The user can open alternative visualization contexts with
48: . PetscViewerASCIIOpen() - output to a specified file
50: Level: beginner
52: .keywords: SLES, view
54: .seealso: PetscViewerASCIIOpen()
55: @*/
56: int SLESView(SLES sles,PetscViewer viewer)
57: {
58: KSP ksp;
59: PC pc;
60: int ierr;
64: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(sles->comm);
67: SLESGetPC(sles,&pc);
68: SLESGetKSP(sles,&ksp);
69: KSPView(ksp,viewer);
70: PCView(pc,viewer);
71: return(0);
72: }
76: /*@C
77: SLESSetOptionsPrefix - Sets the prefix used for searching for all
78: SLES options in the database.
80: Collective on SLES
82: Input Parameter:
83: + sles - the SLES context
84: - prefix - the prefix to prepend to all option names
86: Notes:
87: A hyphen (-) must NOT be given at the beginning of the prefix name.
88: The first character of all runtime options is AUTOMATICALLY the
89: hyphen.
91: This prefix is particularly useful for nested use of SLES. For
92: example, the block Jacobi and block diagonal preconditioners use
93: the prefix "sub_" for options relating to the individual blocks.
95: Level: intermediate
97: .keywords: SLES, set, options, prefix, database
99: .seealso: SLESAppendOptionsPrefix(), SLESGetOptionsPrefix()
100: @*/
101: int SLESSetOptionsPrefix(SLES sles,const char prefix[])
102: {
107: PetscObjectSetOptionsPrefix((PetscObject)sles,prefix);
108: KSPSetOptionsPrefix(sles->ksp,prefix);
109: PCSetOptionsPrefix(sles->pc,prefix);
110: return(0);
111: }
115: /*@C
116: SLESAppendOptionsPrefix - Appends to the prefix used for searching for all
117: SLES options in the database.
119: Collective on SLES
121: Input Parameter:
122: + sles - the SLES context
123: - prefix - the prefix to prepend to all option names
125: Notes:
126: A hyphen (-) must NOT be given at the beginning of the prefix name.
127: The first character of all runtime options is AUTOMATICALLY the
128: hyphen.
130: This prefix is particularly useful for nested use of SLES. For
131: example, the block Jacobi and block diagonal preconditioners use
132: the prefix "sub_" for options relating to the individual blocks.
134: Level: intermediate
136: .keywords: SLES, append, options, prefix, database
138: .seealso: SLESSetOptionsPrefix(), SLESGetOptionsPrefix()
139: @*/
140: int SLESAppendOptionsPrefix(SLES sles,const char prefix[])
141: {
146: PetscObjectAppendOptionsPrefix((PetscObject)sles,prefix);
147: KSPAppendOptionsPrefix(sles->ksp,prefix);
148: PCAppendOptionsPrefix(sles->pc,prefix);
149: return(0);
150: }
154: /*@C
155: SLESGetOptionsPrefix - Gets the prefix used for searching for all
156: SLES options in the database.
158: Not Collective
160: Input Parameter:
161: . sles - the SLES context
163: Output Parameter:
164: . prefix - pointer to the prefix string used
166: Notes:
167: This prefix is particularly useful for nested use of SLES. For
168: example, the block Jacobi and block diagonal preconditioners use
169: the prefix "sub" for options relating to the individual blocks.
171: On the fortran side, the user should pass in a string 'prifix' of
172: sufficient length to hold the prefix.
174: Level: intermediate
176: .keywords: SLES, get, options, prefix, database
178: .seealso: SLESSetOptionsPrefix()
179: @*/
180: int SLESGetOptionsPrefix(SLES sles,char *prefix[])
181: {
186: PetscObjectGetOptionsPrefix((PetscObject)sles,prefix);
187: return(0);
188: }
192: /*@
193: SLESSetFromOptions - Sets various SLES parameters from user options.
194: Also takes all KSP and PC options.
196: Collective on SLES
198: Input Parameter:
199: . sles - the SLES context
201: Level: beginner
203: Options Database Keys:
204: + -sles_view - print actual solver options used
205: - -sles_view_binary - save matrix and vector used for solver
207: .keywords: SLES, set, options, database
209: .seealso: KSPSetFromOptions(),
210: PCSetFromOptions(), SLESGetPC(), SLESGetKSP()
211: @*/
212: int SLESSetFromOptions(SLES sles)
213: {
214: PC tmp_pc;
215: int ierr;
216: PetscTruth flag;
220: PetscOptionsBegin(sles->comm,sles->prefix,"Linear solver (SLES) options","SLES");
221: PetscOptionsName("-sles_diagonal_scale","Diagonal scale matrix before building preconditioner","SLESSetDiagonalScale",&flag);
222: if (flag) {
223: SLESSetDiagonalScale(sles,PETSC_TRUE);
224: }
225: PetscOptionsName("-sles_diagonal_scale_fix","Fix diagonaled scaled matrix after solve","SLESSetDiagonalScaleFix",&flag);
226: if (flag) {
227: SLESSetDiagonalScaleFix(sles);
228: }
229: PetscOptionsName("-sles_view","Print detailed information on solver used","SLESView",0);
230: PetscOptionsName("-sles_view_binary","Saves matrix and vector used in solve","SLESSetFromOptions",0);
231: PetscOptionsEnd();
232: KSPGetPC(sles->ksp,&tmp_pc);
233: if (tmp_pc && sles->pc && tmp_pc!=sles->pc)
234: SETERRQ(1,"Contradictory PCs in SLES and KSP");
235: KSPSetPC(sles->ksp,sles->pc);
236: KSPSetFromOptions(sles->ksp);
237: PCSetFromOptions(sles->pc);
238: return(0);
239: }
243: /*@C
244: SLESCreate - Creates a linear equation solver context.
246: Collective on MPI_Comm
248: Input Parameter:
249: . comm - MPI communicator
251: Output Parameter:
252: . sles - the newly created SLES context
254: Level: beginner
256: .keywords: SLES, create, context
258: .seealso: SLESSolve(), SLESDestroy(), SLES
259: @*/
260: int SLESCreate(MPI_Comm comm,SLES *outsles)
261: {
262: int ierr;
263: SLES sles;
267: *outsles = 0;
268: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
269: SLESInitializePackage(PETSC_NULL);
270: #endif
272: PetscHeaderCreate(sles,_p_SLES,int,SLES_COOKIE,0,"SLES",comm,SLESDestroy,SLESView);
273: PetscLogObjectCreate(sles);
274: sles->bops->publish = SLESPublish_Petsc;
275: KSPCreate(comm,&sles->ksp);
276: PCCreate(comm,&sles->pc);
277: PetscLogObjectParent(sles,sles->ksp);
278: PetscLogObjectParent(sles,sles->pc);
279: sles->setupcalled = 0;
280: sles->dscale = PETSC_FALSE;
281: sles->dscalefix = PETSC_FALSE;
282: sles->diagonal = PETSC_NULL;
283: *outsles = sles;
284: PetscPublishAll(sles);
285: return(0);
286: }
290: /*@C
291: SLESDestroy - Destroys the SLES context.
293: Collective on SLES
295: Input Parameters:
296: . sles - the SLES context
298: Level: beginner
300: .keywords: SLES, destroy, context
302: .seealso: SLESCreate(), SLESSolve()
303: @*/
304: int SLESDestroy(SLES sles)
305: {
310: if (--sles->refct > 0) return(0);
312: PetscObjectDepublish(sles);
313: KSPDestroy(sles->ksp);
314: PCDestroy(sles->pc);
315: if (sles->diagonal) {VecDestroy(sles->diagonal);}
316: PetscLogObjectDestroy(sles);
317: PetscHeaderDestroy(sles);
318: return(0);
319: }
323: /*@
324: SLESSetUp - Performs set up required for solving a linear system.
326: Collective on SLES
328: Input Parameters:
329: + sles - the SLES context
330: . b - the right hand side
331: - x - location to hold solution
333: Note:
334: For basic use of the SLES solvers the user need not explicitly call
335: SLESSetUp(), since these actions will automatically occur during
336: the call to SLESSolve(). However, if one wishes to generate
337: performance data for this computational phase (for example, for
338: incomplete factorization using the ILU preconditioner) using the
339: PETSc log facilities, calling SLESSetUp() is required.
341: Level: advanced
343: .keywords: SLES, solve, linear system
345: .seealso: SLESCreate(), SLESDestroy(), SLESDestroy(), SLESSetUpOnBlocks()
346: @*/
347: int SLESSetUp(SLES sles,Vec b,Vec x)
348: {
350: KSP ksp;
351: PC pc,tmp_pc;
352: Mat mat,pmat;
361: ksp = sles->ksp; pc = sles->pc;
362: KSPSetRhs(ksp,b);
363: KSPSetSolution(ksp,x);
364: if (!sles->setupcalled) {
365: PetscLogEventBegin(SLES_SetUp,sles,b,x,0);
366: KSPGetPC(ksp,&tmp_pc);
367: if (tmp_pc && pc && tmp_pc!=pc)
368: SETERRQ(1,"Contradictory PCs in SLES and KSP");
369: KSPSetPC(ksp,pc);
370: PCSetVector(pc,b);
371: KSPSetUp(sles->ksp);
373: /* scale the matrix if requested */
374: if (sles->dscale) {
375: PCGetOperators(pc,&mat,&pmat,PETSC_NULL);
376: if (mat == pmat) {
377: PetscScalar *xx;
378: int i,n;
379: PetscTruth zeroflag = PETSC_FALSE;
382: if (!sles->diagonal) { /* allocate vector to hold diagonal */
383: VecDuplicate(b,&sles->diagonal);
384: }
385: MatGetDiagonal(mat,sles->diagonal);
386: VecGetLocalSize(sles->diagonal,&n);
387: VecGetArray(sles->diagonal,&xx);
388: for (i=0; i<n; i++) {
389: if (xx[i] != 0.0) xx[i] = 1.0/sqrt(PetscAbsScalar(xx[i]));
390: else {
391: xx[i] = 1.0;
392: zeroflag = PETSC_TRUE;
393: }
394: }
395: VecRestoreArray(sles->diagonal,&xx);
396: if (zeroflag) {
397: PetscLogInfo(pc,"SLESSetUp:Zero detected in diagonal of matrix, using 1 at those locations\n");
398: }
399: MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
400: sles->dscalefix2 = PETSC_FALSE;
401: }
402: }
404: PCSetUp(sles->pc);
405: sles->setupcalled = 1;
406: PetscLogEventEnd(SLES_SetUp,sles,b,x,0);
407: }
408: return(0);
409: }
413: /*@
414: SLESSolve - Solves a linear system.
416: Collective on SLES
418: Input Parameters:
419: + sles - the SLES context
420: - b - the right hand side
422: Output Parameters:
423: + x - the approximate solution
424: - its - the number of iterations until termination
426: Options Database Keys:
427: + -sles_diagonal_scale - diagonally scale linear system before solving
428: . -sles_diagonal_scale_fix - unscale the matrix and right hand side when done
429: . -sles_view - print information about the solver used
430: - -sles_view_binary - save the matrix and right hand side to the default binary file (can be
431: read later with src/sles/examples/tutorials/ex10.c for testing solvers)
433: Notes:
434: Call KSPGetConvergedReason() to determine if the solver converged or failed and
435: why.
437: On return, the parameter "its" contains the iteration number at which convergence
438: was successfully reached or divergence/failure was detected
440: If using a direct method (e.g., via the KSP solver
441: KSPPREONLY and a preconditioner such as PCLU/PCILU),
442: then its=1. See KSPSetTolerances() and KSPDefaultConverged()
443: for more details. Also, see KSPSolve() for more details
444: about iterative solver options.
446: Setting a Nonzero Initial Guess:
447: By default, SLES assumes an initial guess of zero by zeroing
448: the initial value for the solution vector, x. To use a nonzero
449: initial guess, the user must call
450: .vb
451: SLESGetKSP(sles,&ksp);
452: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
453: .ve
455: Solving Successive Linear Systems:
456: When solving multiple linear systems of the same size with the
457: same method, several options are available.
459: (1) To solve successive linear systems having the SAME
460: preconditioner matrix (i.e., the same data structure with exactly
461: the same matrix elements) but different right-hand-side vectors,
462: the user should simply call SLESSolve() multiple times. The
463: preconditioner setup operations (e.g., factorization for ILU) will
464: be done during the first call to SLESSolve() only; such operations
465: will NOT be repeated for successive solves.
467: (2) To solve successive linear systems that have DIFFERENT
468: preconditioner matrices (i.e., the matrix elements and/or the
469: matrix data structure change), the user MUST call SLESSetOperators()
470: and SLESSolve() for each solve. See SLESSetOperators() for
471: options that can save work for such cases.
473: Level: beginner
475: .keywords: SLES, solve, linear system
477: .seealso: SLESCreate(), SLESSetOperators(), SLESGetKSP(), KSPSetTolerances(),
478: KSPDefaultConverged(), KSPSetInitialGuessNonzero(), KSPGetResidualNorm(),
479: KSPSolve()
480: @*/
481: int SLESSolve(SLES sles,Vec b,Vec x,int *its)
482: {
483: int ierr;
484: PetscTruth flg;
485: KSP ksp;
486: PC pc;
494: if (b == x) SETERRQ(PETSC_ERR_ARG_IDN,"b and x must be different vectors");
496: PetscOptionsHasName(sles->prefix,"-sles_view_binary",&flg);
497: if (flg) {
498: Mat mat;
499: PCGetOperators(sles->pc,&mat,PETSC_NULL,PETSC_NULL);
500: MatView(mat,PETSC_VIEWER_BINARY_(sles->comm));
501: VecView(b,PETSC_VIEWER_BINARY_(sles->comm));
502: }
503: ksp = sles->ksp;
504: pc = sles->pc;
505: SLESSetUp(sles,b,x);
506: SLESSetUpOnBlocks(sles);
507: PetscLogEventBegin(SLES_Solve,sles,b,x,0);
508: /* diagonal scale RHS if called for */
509: if (sles->dscale) {
510: VecPointwiseMult(sles->diagonal,b,b);
511: /* second time in, but matrix was scaled back to original */
512: if (sles->dscalefix && sles->dscalefix2) {
513: Mat mat;
515: PCGetOperators(pc,&mat,PETSC_NULL,PETSC_NULL);
516: MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
517: }
518: }
519: PCPreSolve(pc,ksp);
520: KSPSolve(ksp,its);
521: PCPostSolve(pc,ksp);
522: /* diagonal scale solution if called for */
523: if (sles->dscale) {
524: VecPointwiseMult(sles->diagonal,x,x);
525: /* unscale right hand side and matrix */
526: if (sles->dscalefix) {
527: Mat mat;
529: VecReciprocal(sles->diagonal);
530: VecPointwiseMult(sles->diagonal,b,b);
531: PCGetOperators(pc,&mat,PETSC_NULL,PETSC_NULL);
532: MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
533: VecReciprocal(sles->diagonal);
534: sles->dscalefix2 = PETSC_TRUE;
535: }
536: }
537: PetscLogEventEnd(SLES_Solve,sles,b,x,0);
538: PetscOptionsHasName(sles->prefix,"-sles_view",&flg);
539: if (flg && !PetscPreLoadingOn) { SLESView(sles,PETSC_VIEWER_STDOUT_(sles->comm)); }
540: return(0);
541: }
545: /*@
546: SLESSolveTranspose - Solves the transpose of a linear system.
548: Collective on SLES
550: Input Parameters:
551: + sles - the SLES context
552: - b - the right hand side
554: Output Parameters:
555: + x - the approximate solution
556: - its - the number of iterations until termination
558: Notes:
559: On return, the parameter "its" contains
560: + <its> - the iteration number at which convergence
561: was successfully reached, or
562: - <-its> - the negative of the iteration at which
563: divergence or breakdown was detected.
565: Currently works only with the KSPPREONLY method.
567: Setting a Nonzero Initial Guess:
568: By default, SLES assumes an initial guess of zero by zeroing
569: the initial value for the solution vector, x. To use a nonzero
570: initial guess, the user must call
571: .vb
572: SLESGetKSP(sles,&ksp);
573: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
574: .ve
576: Solving Successive Linear Systems:
577: When solving multiple linear systems of the same size with the
578: same method, several options are available.
580: (1) To solve successive linear systems having the SAME
581: preconditioner matrix (i.e., the same data structure with exactly
582: the same matrix elements) but different right-hand-side vectors,
583: the user should simply call SLESSolve() multiple times. The
584: preconditioner setup operations (e.g., factorization for ILU) will
585: be done during the first call to SLESSolve() only; such operations
586: will NOT be repeated for successive solves.
588: (2) To solve successive linear systems that have DIFFERENT
589: preconditioner matrices (i.e., the matrix elements and/or the
590: matrix data structure change), the user MUST call SLESSetOperators()
591: and SLESSolve() for each solve. See SLESSetOperators() for
592: options that can save work for such cases.
594: Level: intermediate
596: .keywords: SLES, solve, linear system
598: .seealso: SLESCreate(), SLESSetOperators(), SLESGetKSP(), KSPSetTolerances(),
599: KSPDefaultConverged(), KSPSetInitialGuessNonzero(), KSPGetResidualNorm(),
600: KSPSolve()
601: @*/
602: int SLESSolveTranspose(SLES sles,Vec b,Vec x,int *its)
603: {
604: int ierr;
605: PetscTruth flg;
606: KSP ksp;
607: PC pc,tmp_pc;
613: if (b == x) SETERRQ(PETSC_ERR_ARG_IDN,"b and x must be different vectors");
617: ksp = sles->ksp; pc = sles->pc;
618: KSPSetRhs(ksp,b);
619: KSPSetSolution(ksp,x);
620: KSPGetPC(ksp,&tmp_pc);
621: if (tmp_pc && pc && tmp_pc!=pc)
622: SETERRQ(1,"Contradictory PCs in SLES and KSP");
623: KSPSetPC(ksp,pc);
624: if (!sles->setupcalled) {
625: SLESSetUp(sles,b,x);
626: }
627: PCSetUpOnBlocks(pc);
628: PetscLogEventBegin(SLES_Solve,sles,b,x,0);
629: KSPSolveTranspose(ksp,its);
630: PetscLogEventEnd(SLES_Solve,sles,b,x,0);
631: PetscOptionsHasName(sles->prefix,"-sles_view",&flg);
632: if (flg) { SLESView(sles,PETSC_VIEWER_STDOUT_(sles->comm)); }
633: return(0);
634: }
638: /*@C
639: SLESGetKSP - Returns the KSP context for a SLES solver.
641: Not Collective, but if KSP will be a parallel object if SLES is
643: Input Parameter:
644: . sles - the SLES context
646: Output Parameter:
647: . ksp - the Krylov space context
649: Notes:
650: The user can then directly manipulate the KSP context to set various
651: options (e.g., by calling KSPSetType()), etc.
653: Level: beginner
654:
655: .keywords: SLES, get, KSP, context
657: .seealso: SLESGetPC()
658: @*/
659: int SLESGetKSP(SLES sles,KSP *ksp)
660: {
663: *ksp = sles->ksp;
664: return(0);
665: }
669: /*@C
670: SLESGetPC - Returns the preconditioner (PC) context for a SLES solver.
672: Not Collective, but if PC will be a parallel object if SLES is
674: Input Parameter:
675: . sles - the SLES context
677: Output Parameter:
678: . pc - the preconditioner context
680: Notes:
681: The user can then directly manipulate the PC context to set various
682: options (e.g., by calling PCSetType()), etc.
684: Level: beginner
686: .keywords: SLES, get, PC, context
688: .seealso: SLESGetKSP()
689: @*/
690: int SLESGetPC(SLES sles,PC *pc)
691: {
694: *pc = sles->pc;
695: return(0);
696: }
698: #include src/mat/matimpl.h
701: /*@
702: SLESSetOperators - Sets the matrix associated with the linear system
703: and a (possibly) different one associated with the preconditioner.
705: Collective on SLES and Mat
707: Input Parameters:
708: + sles - the SLES context
709: . Amat - the matrix associated with the linear system
710: . Pmat - the matrix to be used in constructing the preconditioner, usually the
711: same as Amat.
712: - flag - flag indicating information about the preconditioner matrix structure
713: during successive linear solves. This flag is ignored the first time a
714: linear system is solved, and thus is irrelevant when solving just one linear
715: system.
717: Notes:
718: The flag can be used to eliminate unnecessary work in the preconditioner
719: during the repeated solution of linear systems of the same size. The
720: available options are
721: $ SAME_PRECONDITIONER -
722: $ Pmat is identical during successive linear solves.
723: $ This option is intended for folks who are using
724: $ different Amat and Pmat matrices and want to reuse the
725: $ same preconditioner matrix. For example, this option
726: $ saves work by not recomputing incomplete factorization
727: $ for ILU/ICC preconditioners.
728: $ SAME_NONZERO_PATTERN -
729: $ Pmat has the same nonzero structure during
730: $ successive linear solves.
731: $ DIFFERENT_NONZERO_PATTERN -
732: $ Pmat does not have the same nonzero structure.
734: Caution:
735: If you specify SAME_NONZERO_PATTERN, PETSc believes your assertion
736: and does not check the structure of the matrix. If you erroneously
737: claim that the structure is the same when it actually is not, the new
738: preconditioner will not function correctly. Thus, use this optimization
739: feature carefully!
741: If in doubt about whether your preconditioner matrix has changed
742: structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
744: Level: beginner
746: .keywords: SLES, set, operators, matrix, preconditioner, linear system
748: .seealso: SLESSolve(), SLESGetPC(), PCGetOperators()
749: @*/
750: int SLESSetOperators(SLES sles,Mat Amat,Mat Pmat,MatStructure flag)
751: {
756: PCSetOperators(sles->pc,Amat,Pmat,flag);
757: sles->setupcalled = 0; /* so that next solve call will call setup */
758: return(0);
759: }
763: /*@
764: SLESSetUpOnBlocks - Sets up the preconditioner for each block in
765: the block Jacobi, block Gauss-Seidel, and overlapping Schwarz
766: methods.
768: Collective on SLES
770: Input Parameter:
771: . sles - the SLES context
773: Notes:
774: SLESSetUpOnBlocks() is a routine that the user can optinally call for
775: more precise profiling (via -log_summary) of the setup phase for these
776: block preconditioners. If the user does not call SLESSetUpOnBlocks(),
777: it will automatically be called from within SLESSolve().
778:
779: Calling SLESSetUpOnBlocks() is the same as calling PCSetUpOnBlocks()
780: on the PC context within the SLES context.
782: Level: advanced
784: .keywords: SLES, setup, blocks
786: .seealso: PCSetUpOnBlocks(), SLESSetUp(), PCSetUp()
787: @*/
788: int SLESSetUpOnBlocks(SLES sles)
789: {
791: PC pc;
795: SLESGetPC(sles,&pc);
796: PCSetUpOnBlocks(pc);
797: return(0);
798: }
802: /*@C
803: SLESSetDiagonalScale - Tells SLES to diagonally scale the system
804: before solving. This actually CHANGES the matrix (and right hand side).
806: Collective on SLES
808: Input Parameter:
809: + sles - the SLES context
810: - scale - PETSC_TRUE or PETSC_FALSE
812: Notes:
813: BE CAREFUL with this routine: it actually scales the matrix and right
814: hand side that define the system. After the system is solved the matrix
815: and right hand side remain scaled.
817: This routine is only used if the matrix and preconditioner matrix are
818: the same thing.
819:
820: If you use this with the PCType Eisenstat preconditioner than you can
821: use the PCEisenstatNoDiagonalScaling() option, or -pc_eisenstat_no_diagonal_scaling
822: to save some unneeded, redundant flops.
824: Level: intermediate
826: .keywords: SLES, set, options, prefix, database
828: .seealso: SLESGetDiagonalScale(), SLESSetDiagonalScaleFix()
829: @*/
830: int SLESSetDiagonalScale(SLES sles,PetscTruth scale)
831: {
834: sles->dscale = scale;
835: return(0);
836: }
840: /*@C
841: SLESGetDiagonalScale - Checks if SLES solver scales the matrix and
842: right hand side
844: Not Collective
846: Input Parameter:
847: . sles - the SLES context
849: Output Parameter:
850: . scale - PETSC_TRUE or PETSC_FALSE
852: Notes:
853: BE CAREFUL with this routine: it actually scales the matrix and right
854: hand side that define the system. After the system is solved the matrix
855: and right hand side remain scaled.
857: This routine is only used if the matrix and preconditioner matrix are
858: the same thing.
860: Level: intermediate
862: .keywords: SLES, set, options, prefix, database
864: .seealso: SLESSetDiagonalScale(), SLESSetDiagonalScaleFix()
865: @*/
866: int SLESGetDiagonalScale(SLES sles,PetscTruth *scale)
867: {
870: *scale = sles->dscale;
871: return(0);
872: }
876: /*@C
877: SLESSetDiagonalScaleFix - Tells SLES to diagonally scale the system
878: back after solving.
880: Collective on SLES
882: Input Parameter:
883: . sles - the SLES context
886: Notes:
887: Must be called after SLESSetDiagonalScale()
889: Using this will slow things down, because it rescales the matrix before and
890: after each linear solve. This is intended mainly for testing to allow one
891: to easily get back the original system to make sure the solution computed is
892: accurate enough.
894: This routine is only used if the matrix and preconditioner matrix are
895: the same thing.
897: Level: intermediate
899: .keywords: SLES, set, options, prefix, database
901: .seealso: SLESGetDiagonalScale(), SLESSetDiagonalScale()
902: @*/
903: int SLESSetDiagonalScaleFix(SLES sles)
904: {
907: if (!sles->dscale) {
908: SETERRQ(1,"Must call after SLESSetDiagonalScale(); if setting with command line options you must set\n\
909: -sles_diagonal_scale option if you set the -sles_diagonal_scale_fix option");
910: }
911: sles->dscalefix = PETSC_TRUE;
912: return(0);
913: }