MatGalerkin#
Constructs the coarse grid problem matrix via Galerkin projection.
Synopsis#
#include "petscmat.h"
PetscErrorCode MatGalerkin(Mat restrct, Mat dA, Mat interpolate, MatReuse reuse, PetscReal fill, Mat *A)
If the interpolation and restriction operators are the same, uses MatPtAP()
.
If they are not the same, uses MatMatMatMult()
.
Once the coarse grid problem is constructed, correct for interpolation operators that are not of full rank, which can legitimately happen in the case of non-nested geometric multigrid.
Input Parameters#
restrct - restriction operator
dA - fine grid matrix
interpolate - interpolation operator
reuse - either
MAT_INITIAL_MATRIX
orMAT_REUSE_MATRIX
fill - expected fill, use
PETSC_DETERMINE
orPETSC_DETERMINE
if you do not have a good estimate
Output Parameter#
A - the Galerkin coarse matrix
Options Database Key#
-pc_mg_galerkin <both,pmat,mat,none> - for what matrices the Galerkin process should be used
Note#
The deprecated PETSC_DEFAULT
in fill
also means use the current value
See Also#
Matrices, Mat
, MatPtAP()
, MatMatMatMult()
Level#
developer
Location#
Index of all Mat routines
Table of Contents for all manual pages
Index of all manual pages