How to solve the generalized Lyapunov equation of badly conditioned systems?

Dear community,

I am trying to determine the controllability and observability Gramians of large sparse descriptor systems

using the Generalized lyapunov equations, such as

.

MATLAB can either solve the direct Gramian (using lyap) or its cholesky factor (using gram or lyapchol). The problem I encounter is that, although my system is stable upon evaluation of the eigenvalues (all negative real part), these functions all return the error that the system is unstable and therefore cannot be solved.

An alternative approach is to convert to system to standard statespace, by

which does in fact allow me to solve the Lyapunov equation, but requires an additional inverse. In addition, these structural models often have a badly conditioned E matrix, making the inverse unreliable.

Would you have a suggestion on how to approach this problem most succesfully? I've attached some test that I did to compare different approaches for a set of models. The results, in terms of the norm of the residual, are also included.

Thank you!