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: #undef __FUNCT__  
 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: }

 25: #undef __FUNCT__  
 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: }

 74: #undef __FUNCT__  
 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,char *prefix)
102: {

107:   PetscObjectSetOptionsPrefix((PetscObject)sles,prefix);
108:   KSPSetOptionsPrefix(sles->ksp,prefix);
109:   PCSetOptionsPrefix(sles->pc,prefix);
110:   return(0);
111: }

113: #undef __FUNCT__  
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,char *prefix)
141: {

146:   PetscObjectAppendOptionsPrefix((PetscObject)sles,prefix);
147:   KSPAppendOptionsPrefix(sles->ksp,prefix);
148:   PCAppendOptionsPrefix(sles->pc,prefix);
149:   return(0);
150: }

152: #undef __FUNCT__  
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: }

190: #undef __FUNCT__  
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:   int        ierr;
215:   PetscTruth flag;

219:   PetscOptionsBegin(sles->comm,sles->prefix,"Linear solver (SLES) options","SLES");
220:     PetscOptionsName("-sles_diagonal_scale","Diagonal scale matrix before building preconditioner","SLESSetDiagonalScale",&flag);
221:     if (flag) {
222:       SLESSetDiagonalScale(sles,PETSC_TRUE);
223:     }
224:     PetscOptionsName("-sles_diagonal_scale_fix","Fix diagonaled scaled matrix after solve","SLESSetDiagonalScaleFix",&flag);
225:     if (flag) {
226:       SLESSetDiagonalScaleFix(sles);
227:     }
228:     PetscOptionsName("-sles_view","Print detailed information on solver used","SLESView",0);
229:     PetscOptionsName("-sles_view_binary","Saves matrix and vector used in solve","SLESSetFromOptions",0);
230:   PetscOptionsEnd();
231:   KSPSetPC(sles->ksp,sles->pc);
232:   KSPSetFromOptions(sles->ksp);
233:   PCSetFromOptions(sles->pc);
234:   return(0);
235: }

237: #undef __FUNCT__  
239: /*@C
240:    SLESCreate - Creates a linear equation solver context.

242:    Collective on MPI_Comm

244:    Input Parameter:
245: .  comm - MPI communicator

247:    Output Parameter:
248: .  sles - the newly created SLES context

250:    Level: beginner

252: .keywords: SLES, create, context

254: .seealso: SLESSolve(), SLESDestroy(), SLES
255: @*/
256: int SLESCreate(MPI_Comm comm,SLES *outsles)
257: {
258:   int  ierr;
259:   SLES sles;

263:   *outsles = 0;
264: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
265:   SLESInitializePackage(PETSC_NULL);
266: #endif

268:   PetscHeaderCreate(sles,_p_SLES,int,SLES_COOKIE,0,"SLES",comm,SLESDestroy,SLESView);
269:   PetscLogObjectCreate(sles);
270:   sles->bops->publish = SLESPublish_Petsc;
271:   KSPCreate(comm,&sles->ksp);
272:   PCCreate(comm,&sles->pc);
273:   PetscLogObjectParent(sles,sles->ksp);
274:   PetscLogObjectParent(sles,sles->pc);
275:   sles->setupcalled = 0;
276:   sles->dscale      = PETSC_FALSE;
277:   sles->dscalefix   = PETSC_FALSE;
278:   sles->diagonal    = PETSC_NULL;
279:   *outsles = sles;
280:   PetscPublishAll(sles);
281:   return(0);
282: }

284: #undef __FUNCT__  
286: /*@C
287:    SLESDestroy - Destroys the SLES context.

289:    Collective on SLES

291:    Input Parameters:
292: .  sles - the SLES context

294:    Level: beginner

296: .keywords: SLES, destroy, context

298: .seealso: SLESCreate(), SLESSolve()
299: @*/
300: int SLESDestroy(SLES sles)
301: {

306:   if (--sles->refct > 0) return(0);

308:   PetscObjectDepublish(sles);
309:   KSPDestroy(sles->ksp);
310:   PCDestroy(sles->pc);
311:   if (sles->diagonal) {VecDestroy(sles->diagonal);}
312:   PetscLogObjectDestroy(sles);
313:   PetscHeaderDestroy(sles);
314:   return(0);
315: }

317: #undef __FUNCT__  
319: /*@
320:    SLESSetUp - Performs set up required for solving a linear system.

322:    Collective on SLES

324:    Input Parameters:
325: +  sles - the SLES context
326: .  b - the right hand side
327: -  x - location to hold solution

329:    Note:
330:    For basic use of the SLES solvers the user need not explicitly call
331:    SLESSetUp(), since these actions will automatically occur during
332:    the call to SLESSolve().  However, if one wishes to generate
333:    performance data for this computational phase (for example, for
334:    incomplete factorization using the ILU preconditioner) using the 
335:    PETSc log facilities, calling SLESSetUp() is required.

337:    Level: advanced

339: .keywords: SLES, solve, linear system

341: .seealso: SLESCreate(), SLESDestroy(), SLESDestroy(), SLESSetUpOnBlocks()
342: @*/
343: int SLESSetUp(SLES sles,Vec b,Vec x)
344: {
346:   KSP ksp;
347:   PC  pc;
348:   Mat mat,pmat;


357:   ksp = sles->ksp; pc = sles->pc;
358:   KSPSetRhs(ksp,b);
359:   KSPSetSolution(ksp,x);
360:   if (!sles->setupcalled) {
361:     PetscLogEventBegin(SLES_SetUp,sles,b,x,0);
362:     KSPSetPC(ksp,pc);
363:     PCSetVector(pc,b);
364:     KSPSetUp(sles->ksp);

366:     /* scale the matrix if requested */
367:     if (sles->dscale) {
368:       PCGetOperators(pc,&mat,&pmat,PETSC_NULL);
369:       if (mat == pmat) {
370:         PetscScalar  *xx;
371:         int          i,n;
372:         PetscTruth   zeroflag = PETSC_FALSE;


375:         if (!sles->diagonal) { /* allocate vector to hold diagonal */
376:           VecDuplicate(b,&sles->diagonal);
377:         }
378:         MatGetDiagonal(mat,sles->diagonal);
379:         VecGetLocalSize(sles->diagonal,&n);
380:         VecGetArray(sles->diagonal,&xx);
381:         for (i=0; i<n; i++) {
382:           if (xx[i] != 0.0) xx[i] = 1.0/sqrt(PetscAbsScalar(xx[i]));
383:           else {
384:             xx[i]     = 1.0;
385:             zeroflag  = PETSC_TRUE;
386:           }
387:         }
388:         VecRestoreArray(sles->diagonal,&xx);
389:         if (zeroflag) {
390:           PetscLogInfo(pc,"SLESSetUp:Zero detected in diagonal of matrix, using 1 at those locationsn");
391:         }
392:         MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
393:         sles->dscalefix2 = PETSC_FALSE;
394:       }
395:     }

397:     PCSetUp(sles->pc);
398:     sles->setupcalled = 1;
399:     PetscLogEventEnd(SLES_SetUp,sles,b,x,0);
400:   }
401:   return(0);
402: }

404: #undef __FUNCT__  
406: /*@
407:    SLESSolve - Solves a linear system.

409:    Collective on SLES

411:    Input Parameters:
412: +  sles - the SLES context
413: -  b - the right hand side

415:    Output Parameters:
416: +  x - the approximate solution
417: -  its - the number of iterations until termination

419:    Options Database Keys:
420: +   -sles_diagonal_scale - diagonally scale linear system before solving
421: .   -sles_diagonal_scale_fix - unscale the matrix and right hand side when done
422: .   -sles_view - print information about the solver used
423: -   -sles_view_binary - save the matrix and right hand side to the default binary file (can be
424:        read later with src/sles/examples/tutorials/ex10.c for testing solvers)

426:    Notes:
427:      Call KSPGetConvergedReason() to determine if the solver converged or failed and 
428:      why.

430:      On return, the parameter "its" contains the iteration number at which convergence
431:        was successfully reached or divergence/failure was detected

433:      If using a direct method (e.g., via the KSP solver
434:      KSPPREONLY and a preconditioner such as PCLU/PCILU),
435:      then its=1.  See KSPSetTolerances() and KSPDefaultConverged()
436:      for more details.  Also, see KSPSolve() for more details
437:      about iterative solver options.

439:    Setting a Nonzero Initial Guess:
440:      By default, SLES assumes an initial guess of zero by zeroing
441:      the initial value for the solution vector, x. To use a nonzero 
442:      initial guess, the user must call
443: .vb
444:         SLESGetKSP(sles,&ksp);
445:         KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
446: .ve

448:    Solving Successive Linear Systems:
449:      When solving multiple linear systems of the same size with the
450:      same method, several options are available.

452:      (1) To solve successive linear systems having the SAME
453:      preconditioner matrix (i.e., the same data structure with exactly
454:      the same matrix elements) but different right-hand-side vectors,
455:      the user should simply call SLESSolve() multiple times.  The
456:      preconditioner setup operations (e.g., factorization for ILU) will
457:      be done during the first call to SLESSolve() only; such operations
458:      will NOT be repeated for successive solves.

460:      (2) To solve successive linear systems that have DIFFERENT
461:      preconditioner matrices (i.e., the matrix elements and/or the
462:      matrix data structure change), the user MUST call SLESSetOperators()
463:      and SLESSolve() for each solve.  See SLESSetOperators() for
464:      options that can save work for such cases.

466:    Level: beginner

468: .keywords: SLES, solve, linear system

470: .seealso: SLESCreate(), SLESSetOperators(), SLESGetKSP(), KSPSetTolerances(),
471:           KSPDefaultConverged(), KSPSetInitialGuessNonzero(), KSPGetResidualNorm(),
472:           KSPSolve()
473: @*/
474: int SLESSolve(SLES sles,Vec b,Vec x,int *its)
475: {
476:   int        ierr;
477:   PetscTruth flg;
478:   KSP        ksp;
479:   PC         pc;

487:   if (b == x) SETERRQ(PETSC_ERR_ARG_IDN,"b and x must be different vectors");

489:   PetscOptionsHasName(sles->prefix,"-sles_view_binary",&flg);
490:   if (flg) {
491:     Mat mat;
492:     PCGetOperators(sles->pc,&mat,PETSC_NULL,PETSC_NULL);
493:     MatView(mat,PETSC_VIEWER_BINARY_(sles->comm));
494:     VecView(b,PETSC_VIEWER_BINARY_(sles->comm));
495:   }
496:   ksp  = sles->ksp;
497:   pc   = sles->pc;
498:   SLESSetUp(sles,b,x);
499:   SLESSetUpOnBlocks(sles);
500:   PetscLogEventBegin(SLES_Solve,sles,b,x,0);
501:   /* diagonal scale RHS if called for */
502:   if (sles->dscale) {
503:     VecPointwiseMult(sles->diagonal,b,b);
504:     /* second time in, but matrix was scaled back to original */
505:     if (sles->dscalefix && sles->dscalefix2) {
506:       Mat mat;

508:       PCGetOperators(pc,&mat,PETSC_NULL,PETSC_NULL);
509:       MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
510:     }
511:   }
512:   PCPreSolve(pc,ksp);
513:   KSPSolve(ksp,its);
514:   PCPostSolve(pc,ksp);
515:   /* diagonal scale solution if called for */
516:   if (sles->dscale) {
517:     VecPointwiseMult(sles->diagonal,x,x);
518:     /* unscale right hand side and matrix */
519:     if (sles->dscalefix) {
520:       Mat mat;

522:       VecReciprocal(sles->diagonal);
523:       VecPointwiseMult(sles->diagonal,b,b);
524:       PCGetOperators(pc,&mat,PETSC_NULL,PETSC_NULL);
525:       MatDiagonalScale(mat,sles->diagonal,sles->diagonal);
526:       VecReciprocal(sles->diagonal);
527:       sles->dscalefix2 = PETSC_TRUE;
528:     }
529:   }
530:   PetscLogEventEnd(SLES_Solve,sles,b,x,0);
531:   PetscOptionsHasName(sles->prefix,"-sles_view",&flg);
532:   if (flg && !PetscPreLoadingOn) { SLESView(sles,PETSC_VIEWER_STDOUT_(sles->comm)); }
533:   return(0);
534: }

536: #undef __FUNCT__  
538: /*@
539:    SLESSolveTranspose - Solves the transpose of a linear system.

541:    Collective on SLES

543:    Input Parameters:
544: +  sles - the SLES context
545: -  b - the right hand side

547:    Output Parameters:
548: +  x - the approximate solution
549: -  its - the number of iterations until termination

551:    Notes:
552:      On return, the parameter "its" contains
553: +    <its> - the iteration number at which convergence
554:        was successfully reached, or
555: -    <-its> - the negative of the iteration at which
556:         divergence or breakdown was detected.

558:      Currently works only with the KSPPREONLY method.

560:    Setting a Nonzero Initial Guess:
561:      By default, SLES assumes an initial guess of zero by zeroing
562:      the initial value for the solution vector, x. To use a nonzero 
563:      initial guess, the user must call
564: .vb
565:         SLESGetKSP(sles,&ksp);
566:         KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
567: .ve

569:    Solving Successive Linear Systems:
570:      When solving multiple linear systems of the same size with the
571:      same method, several options are available.

573:      (1) To solve successive linear systems having the SAME
574:      preconditioner matrix (i.e., the same data structure with exactly
575:      the same matrix elements) but different right-hand-side vectors,
576:      the user should simply call SLESSolve() multiple times.  The
577:      preconditioner setup operations (e.g., factorization for ILU) will
578:      be done during the first call to SLESSolve() only; such operations
579:      will NOT be repeated for successive solves.

581:      (2) To solve successive linear systems that have DIFFERENT
582:      preconditioner matrices (i.e., the matrix elements and/or the
583:      matrix data structure change), the user MUST call SLESSetOperators()
584:      and SLESSolve() for each solve.  See SLESSetOperators() for
585:      options that can save work for such cases.

587:    Level: intermediate

589: .keywords: SLES, solve, linear system

591: .seealso: SLESCreate(), SLESSetOperators(), SLESGetKSP(), KSPSetTolerances(),
592:           KSPDefaultConverged(), KSPSetInitialGuessNonzero(), KSPGetResidualNorm(),
593:           KSPSolve()
594: @*/
595: int SLESSolveTranspose(SLES sles,Vec b,Vec x,int *its)
596: {
597:   int        ierr;
598:   PetscTruth flg;
599:   KSP        ksp;
600:   PC         pc;

606:   if (b == x) SETERRQ(PETSC_ERR_ARG_IDN,"b and x must be different vectors");

610:   ksp = sles->ksp; pc = sles->pc;
611:   KSPSetRhs(ksp,b);
612:   KSPSetSolution(ksp,x);
613:   KSPSetPC(ksp,pc);
614:   if (!sles->setupcalled) {
615:     SLESSetUp(sles,b,x);
616:   }
617:   PCSetUpOnBlocks(pc);
618:   PetscLogEventBegin(SLES_Solve,sles,b,x,0);
619:   KSPSolveTranspose(ksp,its);
620:   PetscLogEventEnd(SLES_Solve,sles,b,x,0);
621:   PetscOptionsHasName(sles->prefix,"-sles_view",&flg);
622:   if (flg) { SLESView(sles,PETSC_VIEWER_STDOUT_(sles->comm)); }
623:   return(0);
624: }

626: #undef __FUNCT__  
628: /*@C
629:    SLESGetKSP - Returns the KSP context for a SLES solver.

631:    Not Collective, but if KSP will be a parallel object if SLES is

633:    Input Parameter:
634: .  sles - the SLES context

636:    Output Parameter:
637: .  ksp - the Krylov space context

639:    Notes:  
640:    The user can then directly manipulate the KSP context to set various 
641:    options (e.g., by calling KSPSetType()), etc.

643:    Level: beginner
644:    
645: .keywords: SLES, get, KSP, context

647: .seealso: SLESGetPC()
648: @*/
649: int SLESGetKSP(SLES sles,KSP *ksp)
650: {
653:   *ksp = sles->ksp;
654:   return(0);
655: }

657: #undef __FUNCT__  
659: /*@C
660:    SLESGetPC - Returns the preconditioner (PC) context for a SLES solver.

662:    Not Collective, but if PC will be a parallel object if SLES is

664:    Input Parameter:
665: .  sles - the SLES context

667:    Output Parameter:
668: .  pc - the preconditioner context

670:    Notes:  
671:    The user can then directly manipulate the PC context to set various 
672:    options (e.g., by calling PCSetType()), etc.

674:    Level: beginner

676: .keywords: SLES, get, PC, context

678: .seealso: SLESGetKSP()
679: @*/
680: int SLESGetPC(SLES sles,PC *pc)
681: {
684:   *pc = sles->pc;
685:   return(0);
686: }

688:  #include src/mat/matimpl.h
689: #undef __FUNCT__  
691: /*@
692:    SLESSetOperators - Sets the matrix associated with the linear system
693:    and a (possibly) different one associated with the preconditioner. 

695:    Collective on SLES and Mat

697:    Input Parameters:
698: +  sles - the SLES context
699: .  Amat - the matrix associated with the linear system
700: .  Pmat - the matrix to be used in constructing the preconditioner, usually the
701:           same as Amat. 
702: -  flag - flag indicating information about the preconditioner matrix structure
703:    during successive linear solves.  This flag is ignored the first time a
704:    linear system is solved, and thus is irrelevant when solving just one linear
705:    system.

707:    Notes: 
708:    The flag can be used to eliminate unnecessary work in the preconditioner 
709:    during the repeated solution of linear systems of the same size.  The
710:    available options are
711: $    SAME_PRECONDITIONER -
712: $      Pmat is identical during successive linear solves.
713: $      This option is intended for folks who are using
714: $      different Amat and Pmat matrices and want to reuse the
715: $      same preconditioner matrix.  For example, this option
716: $      saves work by not recomputing incomplete factorization
717: $      for ILU/ICC preconditioners.
718: $    SAME_NONZERO_PATTERN -
719: $      Pmat has the same nonzero structure during
720: $      successive linear solves. 
721: $    DIFFERENT_NONZERO_PATTERN -
722: $      Pmat does not have the same nonzero structure.

724:     Caution:
725:     If you specify SAME_NONZERO_PATTERN, PETSc believes your assertion
726:     and does not check the structure of the matrix.  If you erroneously
727:     claim that the structure is the same when it actually is not, the new
728:     preconditioner will not function correctly.  Thus, use this optimization
729:     feature carefully!

731:     If in doubt about whether your preconditioner matrix has changed
732:     structure or not, use the flag DIFFERENT_NONZERO_PATTERN.

734:     Level: beginner

736: .keywords: SLES, set, operators, matrix, preconditioner, linear system

738: .seealso: SLESSolve(), SLESGetPC(), PCGetOperators()
739: @*/
740: int SLESSetOperators(SLES sles,Mat Amat,Mat Pmat,MatStructure flag)
741: {
746:   PCSetOperators(sles->pc,Amat,Pmat,flag);
747:   sles->setupcalled = 0;  /* so that next solve call will call setup */
748:   return(0);
749: }

751: #undef __FUNCT__  
753: /*@
754:    SLESSetUpOnBlocks - Sets up the preconditioner for each block in
755:    the block Jacobi, block Gauss-Seidel, and overlapping Schwarz 
756:    methods.

758:    Collective on SLES

760:    Input Parameter:
761: .  sles - the SLES context

763:    Notes:
764:    SLESSetUpOnBlocks() is a routine that the user can optinally call for
765:    more precise profiling (via -log_summary) of the setup phase for these
766:    block preconditioners.  If the user does not call SLESSetUpOnBlocks(),
767:    it will automatically be called from within SLESSolve().
768:    
769:    Calling SLESSetUpOnBlocks() is the same as calling PCSetUpOnBlocks()
770:    on the PC context within the SLES context.

772:    Level: advanced

774: .keywords: SLES, setup, blocks

776: .seealso: PCSetUpOnBlocks(), SLESSetUp(), PCSetUp()
777: @*/
778: int SLESSetUpOnBlocks(SLES sles)
779: {
781:   PC  pc;

785:   SLESGetPC(sles,&pc);
786:   PCSetUpOnBlocks(pc);
787:   return(0);
788: }

790: #undef __FUNCT__  
792: /*@C
793:    SLESSetDiagonalScale - Tells SLES to diagonally scale the system
794:      before solving. This actually CHANGES the matrix (and right hand side).

796:    Collective on SLES

798:    Input Parameter:
799: +  sles - the SLES context
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.
809:  
810:     If you use this with the PCType Eisenstat preconditioner than you can 
811:     use the PCEisenstatNoDiagonalScaling() option, or -pc_eisenstat_no_diagonal_scaling
812:     to save some unneeded, redundant flops.

814:    Level: intermediate

816: .keywords: SLES, set, options, prefix, database

818: .seealso: SLESGetDiagonalScale(), SLESSetDiagonalScaleFix()
819: @*/
820: int SLESSetDiagonalScale(SLES sles,PetscTruth scale)
821: {
824:   sles->dscale = scale;
825:   return(0);
826: }

828: #undef __FUNCT__  
830: /*@C
831:    SLESGetDiagonalScale - Checks if SLES solver scales the matrix and
832:                           right hand side

834:    Not Collective

836:    Input Parameter:
837: .  sles - the SLES context

839:    Output Parameter:
840: .  scale - PETSC_TRUE or PETSC_FALSE

842:    Notes:
843:     BE CAREFUL with this routine: it actually scales the matrix and right 
844:     hand side that define the system. After the system is solved the matrix
845:     and right hand side remain scaled.

847:     This routine is only used if the matrix and preconditioner matrix are
848:     the same thing.

850:    Level: intermediate

852: .keywords: SLES, set, options, prefix, database

854: .seealso: SLESSetDiagonalScale(), SLESSetDiagonalScaleFix()
855: @*/
856: int SLESGetDiagonalScale(SLES sles,PetscTruth *scale)
857: {
860:   *scale = sles->dscale;
861:   return(0);
862: }

864: #undef __FUNCT__  
866: /*@C
867:    SLESSetDiagonalScaleFix - Tells SLES to diagonally scale the system
868:      back after solving.

870:    Collective on SLES

872:    Input Parameter:
873: .  sles - the SLES context


876:    Notes:
877:      Must be called after SLESSetDiagonalScale()

879:      Using this will slow things down, because it rescales the matrix before and
880:      after each linear solve. This is intended mainly for testing to allow one
881:      to easily get back the original system to make sure the solution computed is
882:      accurate enough.

884:     This routine is only used if the matrix and preconditioner matrix are
885:     the same thing.

887:    Level: intermediate

889: .keywords: SLES, set, options, prefix, database

891: .seealso: SLESGetDiagonalScale(), SLESSetDiagonalScale()
892: @*/
893: int SLESSetDiagonalScaleFix(SLES sles)
894: {
897:   if (!sles->dscale) {
898:     SETERRQ(1,"Must call after SLESSetDiagonalScale(); if setting with command line options you must setn
899:                -sles_diagonal_scale option if you set the -sles_diagonal_scale_fix option");
900:   }
901:   sles->dscalefix = PETSC_TRUE;
902:   return(0);
903: }