PCFieldSplitSchurGetSubKSP#
Gets the KSP
contexts used inside the Schur complement based PCFIELDSPLIT
Synopsis#
#include "petscpc.h"
PetscErrorCode PCFieldSplitSchurGetSubKSP(PC pc, PetscInt *n, KSP *subksp[])
Collective
Input Parameter#
pc - the preconditioner context
Output Parameters#
n - the number of splits
subksp - the array of
KSP
contexts
Notes#
After PCFieldSplitSchurGetSubKSP()
the array of KSP
s is to be freed by the user with PetscFree()
(not the KSP
just the array that contains them).
You must call PCSetUp()
before calling PCFieldSplitSchurGetSubKSP()
.
If the fieldsplit type is of type PC_COMPOSITE_SCHUR
, it returns (in order)
1 - the
KSP
used for the (1,1) block2 - the
KSP
used for the Schur complement (not the one used for the interior Schur solver)3 - the
KSP
used for the (1,1) block in the upper triangular factor (if different from that of the (1,1) block).
It returns a null array if the fieldsplit is not of type PC_COMPOSITE_SCHUR
; in this case, you should use PCFieldSplitGetSubKSP()
.
Fortran Notes#
You must pass in a KSP
array that is large enough to contain all the local KSP
s.
You can call PCFieldSplitSchurGetSubKSP
(pc,n,PETSC_NULL_KSP
,ierr) to determine how large the
KSP
array must be.
Developer Notes#
There should be a PCFieldSplitRestoreSubKSP()
instead of requiring the user to call PetscFree()
Should the functionality of PCFieldSplitSchurGetSubKSP()
and PCFieldSplitGetSubKSP()
be merged?
The Fortran interface should be modernized to return directly the array of values.
See Also#
Solving Block Matrices with PCFIELDSPLIT, PC
, PCFIELDSPLIT
, PCFieldSplitSetFields()
, PCFieldSplitSetIS()
, PCFieldSplitGetSubKSP()
Level#
advanced
Location#
Examples#
Implementations#
PCFieldSplitSchurGetSubKSP_FieldSplit() in src/ksp/pc/impls/fieldsplit/fieldsplit.c
Index of all PC routines
Table of Contents for all manual pages
Index of all manual pages