SLESSolve

Solves a linear system.

Synopsis

#include "petscsles.h"    
int SLESSolve(SLES sles,Vec b,Vec x,int *its)
Collective on SLES

Input Parameters

sles - the SLES context
b - the right hand side

Output Parameters

x - the approximate solution
its - the number of iterations until termination

Options Database Keys

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

Notes

Call KSPGetConvergedReason() to determine if the solver converged or failed and why.

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

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

Setting a Nonzero Initial Guess

By default, SLES assumes an initial guess of zero by zeroing the initial value for the solution vector, x. To use a nonzero initial guess, the user must call
        SLESGetKSP(sles,&ksp);
        KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);

Solving Successive Linear Systems

When solving multiple linear systems of the same size with the same method, several options are available.

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

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

Keywords

SLES, solve, linear system

See Also

SLESCreate(), SLESSetOperators(), SLESGetKSP(), KSPSetTolerances(),
KSPDefaultConverged(), KSPSetInitialGuessNonzero(), KSPGetResidualNorm(), KSPSolve()

Level:beginner
Location:
src/sles/interface/sles.c
Index of all SLES routines
Table of Contents for all manual pages
Index of all manual pages src/sles/examples/tutorials/ex1.c.html
src/sles/examples/tutorials/ex2.c.html
src/sles/examples/tutorials/ex3.c.html
src/sles/examples/tutorials/ex4.c.html
src/sles/examples/tutorials/ex5.c.html
src/sles/examples/tutorials/ex7.c.html
src/sles/examples/tutorials/ex8.c.html
src/sles/examples/tutorials/ex9.c.html
src/sles/examples/tutorials/ex10.c.html
src/sles/examples/tutorials/ex11.c.html