slepc-3.12.0 2019-09-30
Report Typos and Errors

LMEDenseLyapunovFact

Computes the Cholesky factor of the solution of a dense Lyapunov equation with rank-1 right-hand side.

Synopsis

#include "slepclme.h" 
PetscErrorCode LMEDenseLyapunovChol(LME lme,PetscScalar *H,PetscInt m,PetscInt ldh,PetscScalar *r,PetscScalar *L,PetscInt ldl,PetscReal *res)
Logically Collective on lme

Input Parameters

lme  - linear matrix equation solver context
H  - coefficient matrix
m  - problem size
ldh  - leading dimension of H
r  - right-hand side vector
ldl  - leading dimension of L

Output Parameter

L  - Cholesky factor of the solution
res  - residual norm

Note

The Lyapunov equation has the form H*X + X*H' = -r*r', where H represents the leading mxm submatrix of argument H, and the solution X = L*L'.

H is assumed to be in upper Hessenberg form, with dimensions (m+1)xm. The last row is used to compute the residual norm, assuming H and r come from the projection onto an Arnoldi basis.

See Also

LMESolve()

Location: src/lme/interface/lmedense.c
Index of all LME routines
Table of Contents for all manual pages
Index of all manual pages