petsc-3.7.1 2016-05-15
Report Typos and Errors

PCFieldSplitSetSchurPre

Indicates if the Schur complement is preconditioned by a preconditioner constructed by the A11 matrix. Otherwise no preconditioner is used.

Synopsis

#include "petscpc.h" 
PetscErrorCode PCFieldSplitSetSchurPre(PC pc,PCFieldSplitSchurPreType ptype,Mat pre)
Collective on PC Many br

Input Parameters

pc - the preconditioner context Many br
ptype - which matrix to use for preconditioning the Schur complement: PC_FIELDSPLIT_SCHUR_PRE_A11 (default), PC_FIELDSPLIT_SCHUR_PRE_SELF, PC_FIELDSPLIT_PRE_USER Many br
userpre - matrix to use for preconditioning, or NULL Many br

Options Database

-pc_fieldsplit_schur_precondition <self,selfp,user,a11,full> default is a11 - Many brNotes: Many br
   If ptype is
       user then the preconditioner for the Schur complement is generated from the user provided matrix (pre argument
            to this function).
       a11 then the preconditioner for the Schur complement is generated from the block diagonal part of the preconditioner
            matrix associated with the Schur complement (i.e. A11), not he Schur complement matrix
       full then the preconditioner for the Schur complement is generated from the exact Schur complement matrix representation computed internally by PFIELDSPLIT (this is expensive)
       self the preconditioner for the Schur complement is generated from the symbolic representation of the Schur complement matrix:
            The only preconditioner that currently works with this symbolic respresentation matrix object is the PCLSC
            preconditioner
       selfp then the preconditioning for the Schur complement is generated from an explicitly-assembled approximation Sp = A11 - A10 inv(diag(A00)) A01
            This is only a good preconditioner when diag(A00) is a good preconditioner for A00. Optionally, A00 can be
            lumped before extracting the diagonal: -fieldsplit_1_mat_schur_complement_ainv_type lump

When solving a saddle point problem, where the A11 block is identically zero, using a11 as the ptype only makes sense Many brwith the additional option -fieldsplit_1_pc_type none. Usually for saddle point problems one would use a ptype of self and Many br-fieldsplit_1_pc_type lsc which uses the least squares commutator to compute a preconditioner for the Schur complement. Many br

Many br

See Also

PCFieldSplitGetSchurPre(), PCFieldSplitGetSubKSP(), PCFIELDSPLIT, PCFieldSplitSetFields(), PCFieldSplitSchurPreType,
MatSchurComplementSetAinvType(), PCLSC Many br

Level:intermediate
Location:
src/ksp/pc/impls/fieldsplit/fieldsplit.c
Index of all PC routines
Table of Contents for all manual pages
Index of all manual pages

Examples

src/snes/examples/tutorials/ex70.c.html