MatSchurComplementGetSubMatrices#

Get the individual submatrices in the Schur complement

Synopsis#

#include "petscksp.h" 
PetscErrorCode MatSchurComplementGetSubMatrices(Mat S, Mat *A00, Mat *Ap00, Mat *A01, Mat *A10, Mat *A11)

Collective

Input Parameter#

  • S - matrix obtained with MatCreateSchurComplement() (or equivalent) and implementing the action of A11 - A10 ksp(A00,Ap00) A01

Output Parameters#

  • A00 - the upper-left block of the original matrix A = [A00 A01; A10 A11]

  • Ap00 - preconditioning matrix for use in ksp(A00,Ap00) to approximate the action of A^{-1}

  • A01 - the upper-right block of the original matrix A = [A00 A01; A10 A11]

  • A10 - the lower-left block of the original matrix A = [A00 A01; A10 A11]

  • A11 - (optional) the lower-right block of the original matrix A = [A00 A01; A10 A11]

Note#

A11 is optional, and thus can be NULL.

The reference counts of the submatrices are not increased before they are returned and the matrices should not be modified or destroyed.

See Also#

KSP: Linear System Solvers, MatCreateNormal(), MatMult(), MatCreate(), MatSchurComplementGetKSP(), MatCreateSchurComplement(), MatSchurComplementUpdateSubMatrices()

Level#

intermediate

Location#

src/ksp/ksp/utils/schurm/schurm.c


Edit on GitLab

Index of all KSP routines
Table of Contents for all manual pages
Index of all manual pages